Reliability polynomials for all terminal reliability networks
by Michael Monagan
Abstract:
Routine for the all terminal reliability polynomial for a graph G with probability of edge failure p.
Application Areas/Subjects:
Combinatorics & Graph Theory
Keywords:
reliability polynomial,
>
restart:
Introduction
The main function in this file is the routine relpoly(G,p). It computes the reliability polynomial a(p) for G the so called ``all-terminal reliability network''. Such a network is modelled by a graph G with multiple edges and with a probability p of failure on each edge. The polynomial a(p) gives the probability that the network is spanned in the presence of edge failures, i.e. the probability that each node remains connected to each other node. We give several examples below.
The algorithm used in this routine is a heuristic which works quite well for sparse networks, and appears to be much much (100 times) faster than the algorithm used in the routine networks[spanpoly] in the Maple library.
>
Reliability polynomials
The call relpoly(G,p) computes the reliability polynomial a(p) from the Graph G on N vertices. The input graph G should be input as a list of N polynomials in the vertices X[i], i=1..N, where the c = coeff(G[i],X[j]) gives the multiplicity of the edges from vertex i to j.
The input p is the edge probability and N is the number of vertices.
>
relpoly := proc(G::list(polynom(nonnegint)),p) local n;
n := nops(G);
reliability(G,p,n)
end:
>
reliability := proc(G,p,n) local c,d,k,u,v,w,t,H;
# Choose first vertex of minimum degree
d := n;
for k to nops(G) while d > 0 do
if G[k] <> infinity then
c := nops(indets(G[k]));
if c < d then d := c; u := k fi
fi;
od;
if d = 0 then if n = 1 then RETURN(1) else RETURN(0) fi fi;
v := op(1,op(1,indets(G[u]))); c := coeff(G[u],X[v],1);
if d = 1 then # Isolated node contraction
H := subsop( u=infinity, v=G[v]-c*X[u], G );
RETURN( expand( (1-(1-p)^c) * reliability(H,p,n-1) ) )
fi;
# Delete the edge of multiplicity c between vertices u & v
# (All edges between u and v have failed with probability (1-p)^c)
H := subsop( u=G[u]-c*X[v], v=G[v]-c*X[u], G );
t := (1-p)^c * reliability(H,p,n);
# Contract vertices u & v (at least one edge has not failed)
H := subsop( u=infinity, v=H[u]+H[v], subs(X[u]=X[v],H) );
t := expand( (1-(1-p)^c) * reliability(H,p,n-1) + t );
# Now remember this computation
reliability(G,p,n) := t;
end:
This is a utility routine to insert the edge u v into the graph G.
G should be initialized to vector (array) of N zeroes.
The third argument m specifies the multiplicity of the edge
insert := proc(G,u,v,m)
if nargs = 3 then insert(G,u,v,1) else
G[u] := G[u]+m*X[v];
G[v] := G[v]+m*X[u];
NULL
fi
end:
Example 1: computing the reliability polynomials for cliques
relpoly_clique(N,p); computes the reliability polynomial for a clique on N vertices with edge failure probability p
>
relpoly_clique := proc(N,p) local G,i,j;
G := array([0$N]);
for i to N-1 do for j from i+1 to N do insert(G,i,j) od od;
G := convert(G,list):
relpoly(G,p);
end:
>
for i from 1 to 6 do
printf(` Clique on %d vertices is\n`,i);
relpoly_clique(i,p);
od;
Clique on 1 vertices is
Clique on 2 vertices is
Clique on 3 vertices is
Clique on 4 vertices is
Clique on 5 vertices is
Clique on 6 vertices is
Example 2: reliability polynomials for simple loops
relpoly_loop(N,p); computes the reliability polynomial for
a simply cycle on N vertices with edge failure probability p
>
relpoly_loop := proc(N,p) local G,i;
G := array([0$N]);
for i to N-1 do insert(G,i,i+1) od;
insert(G,N,1);
G := convert(G,list):
relpoly(G,p);
end:
>
for i from 1 to 6 do
printf( ` Loop of %d vertices\n`, i );
relpoly_clique(i,p);
od;
Loop of 1 vertices
Loop of 2 vertices
Loop of 3 vertices
Loop of 4 vertices
Loop of 5 vertices
Loop of 6 vertices
Example 3: reliability polynomial for the Red Arpanet network
>
relpoly_RedArpanet := proc(p) local N,G;
N := 20;
G := array([0$N]):
insert(G,1,2):
insert(G,1,3):
insert(G,2,3):
insert(G,2,4):
insert(G,3,5):
insert(G,4,5):
insert(G,4,6):
insert(G,6,7):
insert(G,5,7):
insert(G,1,8):
insert(G,8,9):
insert(G,9,10):
insert(G,9,11):
insert(G,10,11):
insert(G,7,10):
insert(G,6,11):
insert(G,8,12):
insert(G,10,15):
insert(G,11,20):
insert(G,12,13):
insert(G,13,14):
insert(G,14,15):
insert(G,13,15):
insert(G,14,18):
insert(G,12,16):
insert(G,15,17):
insert(G,16,17):
insert(G,16,18):
insert(G,17,19):
insert(G,18,19):
insert(G,18,20):
insert(G,19,20):
G := convert(G,list);
relpoly(G,p)
end:
>
relpoly_RedArpanet(p);
>