123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681(**************************************************************************)(* *)(* Ocamlgraph: a generic graph library for OCaml *)(* Copyright (C) 2004-2010 *)(* Sylvain Conchon, Jean-Christophe Filliatre and Julien Signoles *)(* *)(* This software is free software; you can redistribute it and/or *)(* modify it under the terms of the GNU Library General Public *)(* License version 2.1, with the special exception on linking *)(* described in file LICENSE. *)(* *)(* This software is distributed in the hope that it will be useful, *)(* but WITHOUT ANY WARRANTY; without even the implied warranty of *)(* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. *)(* *)(**************************************************************************)(* This module is a contribution of Yuto Takei *)moduleImperative(G:Sig.IM)(W:Sig.WEIGHTwithtypeedge=G.E.t)=structmoduleS=Set.Make(G.V)moduleM=Map.Make(G.V)moduleV=G.VmoduleE=G.E(* [G.t] represents graph itself. [unit M.t] maintains a list of
source vertices to keep track of distances for all vertices.
[(G.E.t option * W.t) M.t M.t] holds mappings for all vertices,
each of which contains its shortest-path tree ancestor (parent)
and a distances from source vertices. *)typet=G.t*S.tref*(G.E.toption*W.t)M.tM.treftypeedge=G.edgetypevertex=G.vertexletsovv=string_of_int(Obj.magic(V.labelv))letdump_cyclecycle=letv0=G.E.src(List.hdcycle)inprint_string("("^(sovv0)^")");letv1=List.fold_left(funve->assert((G.V.comparev(G.E.srce))=0);letv=G.E.dsteinprint_string("-("^(sovv)^")");v)v0cycleinassert(G.V.equalv0v1);print_string"\n"letdump_set=S.iter(funx->print_string((sovx)^", "))letdump(src,dist)=print_string"====================\nS: ";dump_set!src;print_string"\nMap:";M.iter(funkv->print_string("\n "^(sovk)^": ");M.iter(funk(origin,dist)->print_string("("^(sovk)^">>"^(matchoriginwith|None->"---"|Somee->(sov(G.E.srce))^">"^(sov(G.E.dste)))^":"^(string_of_int(Obj.magicdist))^") "))v)!dist;print_string"\n"(* If an edge is going to be added to the graph, which will cause
a negative cycle, raises [Negative_cycle] with edges that can
form such the cycle. *)exceptionNegative_cycleofG.E.tlistletcreate?size()=letg=matchsizewith|Somesize->G.create~size()|None->G.create()in(g,refS.empty,refM.empty)letcopy(g,src,dist)=(G.copyg,ref(!src),ref(!dist))letclear(g,src,dist)=G.clearg;src:=S.empty;dist:=M.emptyletadd_vertex(g,src,dist)v=(* Before adding vertex to the graph, make sure that the vertex
is not in the graph. If already in the graph, just do
nothing and return as is. *)ifnot(G.mem_vertexgv)thenbegin(* Add a vertex to the original one *)G.add_vertexgv;(* The new vertex will immediately be added to the source list *)src:=S.addv!src;(* The new edge should contain a distance mapping with only
from myself with distance zero. *)dist:=M.addv(M.addv(None,W.zero)M.empty)!dist;dump(src,dist)endletrecpropagate(g,src,dist)qstart=ifQueue.is_emptyqthen(g,src,dist)elsebeginlet(v1,v1src)=Queue.popqinletv1dist=M.findv1distinletdist=G.fold_succ_e(funedist->letv2=G.E.dsteinletv2dist=ifM.memv2distthenM.findv2distelseM.emptyin(* Compare distances from given source vertices.
If relax happens, record it to the new list. *)let(v2dist,nextSrc)=S.fold(funx(v2dist,nextSrc)->let_,dev1=M.findxv1distinletndev2=W.adddev1(W.weighte)inletimprovement=trylet_,dev2=M.findxv2distinW.comparendev2dev2<0withNot_found->trueinifimprovementthenletv2dist=M.addx(Somee,ndev2)v2distinletnextSrc=S.addxnextSrcin(v2dist,nextSrc)else(v2dist,nextSrc))v1src(v2dist,S.empty)inifS.is_emptynextSrcthendistelseifG.V.equalstartv2then(* Propagation reaches back to the starting node, which
immediately means presence of a negative cycle. *)(* We should use one of 'src' to traverse to the start node *)letdist=M.addv2v2distdistinletcycle=S.fold(funsx->letrecbuild_cyclexret=matchM.finds(M.findxdist)with|Somee,_->lety=G.E.srceinletcycle=e::retinifG.V.equalstartythenSomecycleelsebuild_cycleycycle|_->Noneinmatchxwith|None->build_cyclev2[]|Some_->x)nextSrcNoneinletcycle=matchcyclewith|Somex->x|None->assertfalseindump_cyclecycle;raise(Negative_cyclecycle)elsebegin(* TODO: Some room for improvement.
If queue has (v2, s) already, technically we can merge
nextSrc into s, so that the number of propagation can be
reduced. *)Queue.push(v2,nextSrc)q;M.addv2v2distdistend)gv1distinpropagate(g,src,dist)qstartendletm_cardinalm=M.fold(fun__acc->acc+1)m0letset_of_mapm=M.fold(funk_acc->S.addkacc)mS.emptyletadd_edge_internal(g,src,dist)v1v2=(* Distance mappings at v1 *)letdv1=M.findv1distin(* To reduce the amount of codes, we just start propagation from
v1. Of course, this can be optimized by starting from v2. But
it may duplicate the same code in multiple places in the
file. In addition, such an optimization only cost for small
amount, which precisely is the operations to relax edges from
v1, other than which have been existed before this
[add_edge_e] call. *)letq=Queue.create()in(* We need to check whether v2 should be kept in the source list
or not. That is, if there maybe a cycle with v1, the
distance from v1 should be still maintained. Otherwise,
simply ignore the distance from v2 *)ifm_cardinaldv1=1&&M.memv2dv1then((* Now we definitely introduced a loop (and possibly non-negative)!
Let me see if this would be negative or not... *)Queue.add(v1,(S.addv2S.empty))q;propagate(g,src,dist)qv1)else((* Or even if we fall back to else-clause here, the edge addition
may have introduced a cycle. Anyway, we need to check if one is
newly created or not at [propagate] *)let(src,dist,dv1)=ifnot(S.memv2src)then(* If v2 isn't one of the source vertices, just simply do
propagation. *)(src,dist,dv1)else(* We can exclude v2 from the list of source because
one can reach v2 from some other vertex. *)letsrc=S.removev2srcin(* Note that following line can be skipped only if the
user don't remove vertex. Otherwise, such operation
like [add_edge g v1 v2] > [remove_vertex g v2] >
[add_vertex g v2] can result in unexpected
behavior. *)letdist=M.map(M.removev2)distin(* We need to re-obtain the distance mappings at v1,
since it can be changed by the line above. *)letdv1=M.findv1distin(src,dist,dv1)in(* Now let's start propagation. *)Queue.add(v1,set_of_mapdv1)q;propagate(g,src,dist)qv1)letadd_edge_e(g,src,dist)e=(* Before adding edge to the graph, make sure that the edge is
not in the graph. If already in the graph, just do nothing
and return as is. *)ifnot(G.mem_edge_ege)thenbegin(* Vertices involved *)letv1=G.E.srceinletv2=G.E.dsteinList.iter(add_vertex(g,src,dist))[v1;v2];begintry(* Because we can restore the graph by calling [G.remove_edge_e]
even in case of failure, we first add it by [G.add_edge_e]. *)G.add_edge_ege;let(_,src',dist')=add_edge_internal(g,!src,!dist)v1v2insrc:=src';dist:=dist'withexp->(* In case of excecption, restore the graph by removing the
edge, and rethrow the exception. *)G.remove_edge_ege;raiseexpend;dump(src,dist)endletadd_edge(g,src,dist)v1v2=(* Same as [add_edge_e] *)ifnot(G.mem_edgegv1v2)thenbeginList.iter(add_vertex(g,src,dist))[v1;v2];begintry(* Because we cannot know the default value for edge length,
we first try to add one by [G.add_edge]. If there occurs an
exception, restore the graph by [G.remove_edge] since there
were no other connections between [v1] and [v2]. *)G.add_edgegv1v2;let(_,src',dist')=add_edge_internal(g,!src,!dist)v1v2insrc:=src';dist:=dist'withexp->(* In case of excecption, restore the graph by removing the
edge, and rethrow the exception. *)G.remove_edgegv1v2;raiseexpend;dump(src,dist)endletremove_edge_internal(g,src)v2=(* Actually, we need to rebuild the distance table, rather than
traverse precedants to remove the edge. *)letq=Queue.create()inprint_string("dump: ");dump_setsrc;letdist=S.fold(funxdist->print_string("source: "^(sovx)^"\n");Queue.add(x,(S.addxS.empty))q;M.addx(M.addx(None,W.zero)M.empty)dist)srcM.emptyinletg,src,dist=propagate(g,src,dist)q(S.choosesrc)inifM.memv2distthen(g,src,dist)else(Queue.add(v2,(S.addv2S.empty))q;letsrc=S.addv2srcinletdist=M.addv2(M.addv2(None,W.zero)M.empty)distinpropagate(g,src,dist)qv2)letremove_edge_e(g,src,dist)e=(* Same as [add_edge_e] *)ifG.mem_edge_egethenbeginG.remove_edge_ege;(* Vertices involved *)letv2=G.E.dsteinlet(_,src',dist')=remove_edge_internal(g,!src)v2insrc:=src';dist:=dist';dump(src,dist)endletremove_edge(g,src,dist)v1v2=(* Same as [add_edge] *)ifG.mem_edgegv1v2thenbeginG.remove_edgegv1v2;let(_,src',dist')=remove_edge_internal(g,!src)v2insrc:=src';dist:=dist';dump(src,dist)endletremove_vertex(g,src,dist)v=(* Same as [add_edge] *)ifG.mem_vertexgvthenbegin(* [remove_vertex] first deletes all outgoing edges from [v] *)G.iter_succ_e(fune->remove_edge_e(g,src,dist)e)gv;(* Then after, deletes all incoming edges to [v] *)G.iter_pred_e(fune->remove_edge_e(g,src,dist)e)gv;(* Note that we are iterating on [g] that is being modified during
iteration. We can do such an above iteration since G is here
permanent. Do not try this for imperative graph. *)(* Now we can feel free to delete [v]. *)G.remove_vertexgv;src:=S.removev!src;dist:=M.removev(M.map(M.removev)!dist);dump(src,dist)endletmap_vertexf(g,src,dist)=letmap_mapupdatem=M.fold(funvmacc->M.add(fv)(updatem)acc)mM.emptyinlet(g,src,dist)=(G.map_vertexfg,S.fold(funvacc->S.add(fv)acc)!srcS.empty,letupdate=function|None,_asv->v|Somee,w->Some(E.create(f(E.srce))(E.labele)(f(E.dste))),winmap_map(map_mapupdate)!dist)in(g,refsrc,refdist)letfold_pred_ef(g,_,_)=G.fold_pred_efgletiter_pred_ef(g,_,_)=G.iter_pred_efgletfold_succ_ef(g,_,_)=G.fold_succ_efgletiter_succ_ef(g,_,_)=G.iter_succ_efgletfold_predf(g,_,_)=G.fold_predfgletfold_succf(g,_,_)=G.fold_succfgletiter_predf(g,_,_)=G.iter_predfgletiter_succf(g,_,_)=G.iter_succfgletfold_edges_ef(g,_,_)=G.fold_edges_efgletiter_edges_ef(g,_,_)=G.iter_edges_efgletfold_edgesf(g,_,_)=G.fold_edgesfgletiter_edgesf(g,_,_)=G.iter_edgesfgletfold_vertexf(g,_,_)=G.fold_vertexfgletiter_vertexf(g,_,_)=G.iter_vertexfgletpred_e(g,_,_)=G.pred_egletsucc_e(g,_,_)=G.succ_egletpred(g,_,_)=G.predgletsucc(g,_,_)=G.succgletfind_all_edges(g,_,_)=G.find_all_edgesgletfind_edge(g,_,_)=G.find_edgegletmem_edge_e(g,_,_)=G.mem_edge_egletmem_edge(g,_,_)=G.mem_edgegletmem_vertex(g,_,_)=G.mem_vertexgletin_degree(g,_,_)=G.in_degreegletout_degree(g,_,_)=G.out_degreegletnb_edges(g,_,_)=G.nb_edgesgletnb_vertex(g,_,_)=G.nb_vertexgletis_empty(g,_,_)=G.is_emptygletis_directed=G.is_directedmoduleMark=structtypegraph=ttypevertex=G.vertexletclearg=let(g,_,_)=ginG.Mark.cleargletget=G.Mark.getletset=G.Mark.setendendmodulePersistent(G:Sig.P)(W:Sig.WEIGHTwithtypeedge=G.E.t)=structmoduleS=Set.Make(G.V)moduleM=Map.Make(G.V)moduleE=G.EmoduleV=G.V(* [G.t] represents graph itself. [unit M.t] maintains a list of
source vertices to keep track of distances for all vertices.
[(G.E.t option * W.t) M.t M.t] holds mappings for all vertices,
each of which contains its shortest-path tree ancestor (parent)
and a distances from source vertices. *)typet=G.t*S.t*(G.E.toption*W.t)M.tM.ttypeedge=G.edgetypevertex=G.vertex(* If an edge is going to be added to the graph, which will cause
a negative cycle, raises [Negative_cycle] with edges that can
form such the cycle. *)exceptionNegative_cycleofG.E.tlistletempty:t=letg=G.emptyinletsrc=S.emptyinletdist=M.emptyin(g,src,dist)letadd_vertex(g,src,dist)v=(* Before adding vertex to the graph, make sure that the vertex
is not in the graph. If already in the graph, just do
nothing and return as is. *)ifG.mem_vertexgvthen(g,src,dist)else(* Add a vertex to the original one *)(G.add_vertexgv),(* The new vertex will immediately be added to the source list *)(S.addvsrc),(* The new edge should contain a distance mapping with only
from myself with distance zero. *)(M.addv(M.addv(None,W.zero)M.empty)dist)letrecpropagate(g,src,dist)qstart=ifQueue.is_emptyqthen(g,src,dist)elsebeginlet(v1,v1src)=Queue.popqinletv1dist=M.findv1distinletdist=G.fold_succ_e(funedist->letv2=G.E.dsteinletv2dist=M.findv2distin(* Compare distances from given source vertices.
If relax happens, record it to the new list. *)let(v2dist,nextSrc)=S.fold(funx(v2dist,nextSrc)->let_,dev1=M.findxv1distinletndev2=W.adddev1(W.weighte)inletimprovement=trylet_,dev2=M.findxv2distinW.comparendev2dev2<0withNot_found->trueinifimprovementthenletv2dist=M.addx(Somee,ndev2)v2distinletnextSrc=S.addxnextSrcin(v2dist,nextSrc)else(v2dist,nextSrc))v1src(v2dist,S.empty)inifS.is_emptynextSrcthendistelseifG.V.equalstartv2then(* Propagation reaches back to the starting node, which
immediately means presence of a negative cycle. *)(* We should use one of 'src' to traverse to the start node *)lets=S.choosenextSrcinletrecbuild_cyclexret=matchM.finds(M.findxdist)with|Somee,_->lety=G.E.srceinletcycle=e::retinifG.V.equalstartythencycleelsebuild_cycleycycle|_->assertfalseinraise(Negative_cycle(build_cyclev2[]))elsebegin(* TODO: Some room for improvement.
If queue has (v2, s) already, technically we can merge
nextSrc into s, so that the number of propagation can be
reduced. *)Queue.push(v2,nextSrc)q;M.addv2v2distdistend)gv1distinpropagate(g,src,dist)qstartendletm_cardinalm=M.fold(fun__acc->acc+1)m0letset_of_mapm=M.fold(funk_acc->S.addkacc)mS.emptyletadd_edge_internal(g,src,dist)v1v2=(* Distance mappings at v1 *)letdv1=M.findv1distin(* To reduce the amount of codes, we just start propagation from
v1. Of course, this can be optimized by starting from v2. But
it may duplicate the same code in multiple places in the
file. In addition, such an optimization only cost for small
amount, which precisely is the operations to relax edges from
v1, other than which have been existed before this
[add_edge_e] call. *)letq=Queue.create()in(* We need to check whether v2 should be kept in the source list
or not. That is, if there maybe a cycle with v1, the
distance from v1 should be still maintained. Otherwise,
simply ignore the distance from v2 *)ifm_cardinaldv1=1&&M.memv2dv1then((* Now we definitely introduced a loop (but possibly non-negative)!
Let me see if this would be negative or not... *)Queue.add(v1,(S.addv2S.empty))q;propagate(g,src,dist)qv1)else((* Or even if we fall back to else-clause here, the edge addition
may have introduced a cycle. Anyway, we need to check if one is
newly created or not at [propagate] *)let(src,dist,dv1)=ifnot(S.memv2src)then(* If v2 isn't one of the source vertices, just simply do
propagation. *)(src,dist,dv1)else(* We can exclude v2 from the list of source because
one can reach v2 from some other vertex. *)((S.removev2src),(* Note that following line can be skipped only if the
user don't remove vertex. Otherwise, such operation
like [add_edge g v1 v2] > [remove_vertex g v2] >
[add_vertex g v2] can result in unexpected
behaviour. *)(M.map(M.removev2)dist),(* We need to re-obtain the distance mappings at v1,
since it can be changed by the line above. *)(M.findv1dist))in(* Now let's start propagation. *)Queue.add(v1,set_of_mapdv1)q;propagate(g,src,dist)qv1)letadd_edge_e(g,src,dist)e=(* Before adding edge to the graph, make sure that the edge is
not in the graph. If already in the graph, just do nothing
and return as is. *)ifG.mem_edge_egethen(g,src,dist)elsebegin(* Vertices involved *)letv1=G.E.srceinletv2=G.E.dsteinlet(g,src,dist)=List.fold_leftadd_vertex(g,src,dist)[v1;v2]inletg=G.add_edge_egeinadd_edge_internal(g,src,dist)v1v2endletadd_edge(g,src,dist)v1v2=(* Same as [add_edge_e] *)ifG.mem_edgegv1v2then(g,src,dist)elsebeginlet(g,src,dist)=List.fold_leftadd_vertex(g,src,dist)[v1;v2]inletg=G.add_edgegv1v2inadd_edge_internal(g,src,dist)v1v2endletremove_edge_internal(g,src)v2=(* Actually, we need to rebuild the distance table, rather than
traverse precedants to remove the edge. *)letq=Queue.create()inletdist=S.fold(funxdist->Queue.add(x,(S.addxS.empty))q;M.addx(M.addx(None,W.zero)M.empty)dist)srcM.emptyinletg,src,dist=propagate(g,src,dist)q(S.choosesrc)inifM.memv2distthen(g,src,dist)else(Queue.add(v2,(S.addv2S.empty))q;letsrc=S.addv2srcinletdist=M.addv2(M.addv2(None,W.zero)M.empty)distinpropagate(g,src,dist)qv2)letremove_edge_e(g,src,dist)e=(* Same as [add_edge_e] *)ifnot(G.mem_edge_ege)then(g,src,dist)elsebeginletg=G.remove_edge_egein(* Vertices involved *)letv2=G.E.dsteinremove_edge_internal(g,src)v2endletremove_edge(g,src,dist)v1v2=(* Same as [add_edge] *)ifnot(G.mem_edgegv1v2)then(g,src,dist)elsebeginletg=G.remove_edgegv1v2inremove_edge_internal(g,src)v2endletremove_vertextv=(* [remove_vertex] first deletes all outgoing edges from [v] *)let(g,_,_)=tinlett=G.fold_succ_e(funet->remove_edge_ete)gvtin(* Then after, deletes all incoming edges to [v] *)let(g,_,_)=tinlett=G.fold_pred_e(funet->remove_edge_ete)gvtin(* Note that we are iterating on [g] that is being modified during
iteration. We can do such an above iteration since G is here
permanent. Do not try this for imperative graph. *)let(g,src,dist)=tin(* Now we can feel free to delete [v]. *)(G.remove_vertexgv,(S.removevsrc),(M.map(M.removev)dist))letmap_vertexf(g,src,dist)=letmap_mapupdatem=M.fold(funvmacc->M.add(fv)(updatem)acc)mM.emptyin(G.map_vertexfg,S.fold(funvacc->S.add(fv)acc)srcS.empty,letupdate=function|None,_asv->v|Somee,w->Some(E.create(f(E.srce))(E.labele)(f(E.dste))),winmap_map(map_mapupdate)dist)(* All below are wrappers *)letfold_pred_ef(g,_,_)=G.fold_pred_efgletiter_pred_ef(g,_,_)=G.iter_pred_efgletfold_succ_ef(g,_,_)=G.fold_succ_efgletiter_succ_ef(g,_,_)=G.iter_succ_efgletfold_predf(g,_,_)=G.fold_predfgletfold_succf(g,_,_)=G.fold_succfgletiter_predf(g,_,_)=G.iter_predfgletiter_succf(g,_,_)=G.iter_succfgletfold_edges_ef(g,_,_)=G.fold_edges_efgletiter_edges_ef(g,_,_)=G.iter_edges_efgletfold_edgesf(g,_,_)=G.fold_edgesfgletiter_edgesf(g,_,_)=G.iter_edgesfgletfold_vertexf(g,_,_)=G.fold_vertexfgletiter_vertexf(g,_,_)=G.iter_vertexfgletpred_e(g,_,_)=G.pred_egletsucc_e(g,_,_)=G.succ_egletpred(g,_,_)=G.predgletsucc(g,_,_)=G.succgletfind_all_edges(g,_,_)=G.find_all_edgesgletfind_edge(g,_,_)=G.find_edgegletmem_edge_e(g,_,_)=G.mem_edge_egletmem_edge(g,_,_)=G.mem_edgegletmem_vertex(g,_,_)=G.mem_vertexgletin_degree(g,_,_)=G.in_degreegletout_degree(g,_,_)=G.out_degreegletnb_edges(g,_,_)=G.nb_edgesgletnb_vertex(g,_,_)=G.nb_vertexgletis_empty(g,_,_)=G.is_emptygletis_directed=G.is_directedend