bvr.m

best routing protocol

function [routing_path, routing_dist, routing_load,storage_load, ...
	  dht_packets, failure_rate, state_per_node, transmission_dist, routing_hops,...
	  nodes_with_2hop_info, extra_control_traffic_2h] = ...
    bvr(neighbor_dist, r, k, beacons, test_pairs, USE_STANDARD_BVR, ...
	failure_nodes, study_dht, run_scoped_flooding, fetch_2hop_info) 
% r is the total number of beacons
% 0 < neighbor_dist(i,j) < Inf: i and j are neighbors, otherwise not
% k: used in distance calculation, k<=r
% beacons: the IDs of chosen beacons. if not given, r nodes are
% randomly picked
% test_pairs: p x 2 matrix p(:,1) is src, p(:,2) is dest

% routing_path: the NxN or p x 1 cell of routing paths
% routing_hops: the NxN or p x 1 matrix of routing hops, Inf means routing
% failure
% routing_load(i) for all pair traffic, the load on vertex i
% fetch_2hop_info: whether to fetch 2-hop neighbor's info on demand
debug = 0;
N = length(neighbor_dist);
if nargin<3
  k = r;
end

if nargin<5
  test_pairs = [];
  p = 0;
else
  [p tmp]= size(test_pairs);
  if tmp~=2
    error('wrong argument with test_pairs');
  end
end

if nargin<6
  USE_STANDARD_BVR = 1;
end

if nargin<7
  failure_nodes = [];
end
if nargin<8
  study_dht = 0;
end

if nargin<9
  run_scoped_flooding = 1; %simulate this
end

if nargin<10
  fetch_2hop_info = 0;
end
tmp = failure_nodes;
failure_nodes = zeros(1,N);
failure_nodes(tmp) = ones(1,length(tmp));
clear tmp;

fprintf('BVR: starting.\n');  
if USE_STANDARD_BVR
  fprintf('BVR: standard distance function\n');
else
  fprintf('BVR: sqr L2 norm\n');
end

if ~issparse(neighbor_dist)
  neighbor_dist = dist_full2sparse(neighbor_dist);
end

nodes_with_2hop_info = zeros(N, 1);
state_per_node = zeros(N,1);
extra_control_traffic_2h = zeros(N,1);
%get neighbor_hops matrix coz that's what is needed in BVR
%connected_mask = (neighbor_dist>0) & (neighbor_dist<Inf)
%neighbor_hops = neighbor_dist - (neighbor_dist-1).*connected_mask
if 0 % make it a little faster as we only deal with hop matrix
  for i=1:N
    for j=1:N
      if i==j
	neighbor_hops(i,j)=0;
      elseif neighbor_dist(i,j)<Inf
	neighbor_hops(i,j)=1;
      else
	neighbor_hops(i,j)=Inf;
      end
    end
  end
else
  neighbor_hops=neighbor_dist;
end

%neighbor_map = neighbor_hops + diag(repmat(Inf,1,N));
neighbor_cell = cell(N,1);
for i=1:N
  neighbor_cell{i}=find(neighbor_hops(i,:));
end

%first pick r beacons randomly if not given
if nargin<4
  tmp = randperm(N);
  beacons = tmp(1:r);
end
beacons = sort(beacons);

%get pairwise shortest path in hops
%[shortest_hops routing_nexthop]= floyd_warshall(neighbor_hops);
%[shortest_hops routing_nexthop]= allpair_dijkstra(neighbor_hops);
%beacon_vectors = shortest_hops(:, beacons);
%node_parents(i,k) is the parent of node i to beacon k
%node_parents = routing_nexthop(:, beacons);
[beacon_vectors node_parents] = allpair_dijkstra(neighbor_hops, beacons);

%C_k(d,:) is the idx of k closest beacons to d, as needed in
%distance calculation
% beacons(idx) is the real index as in neighbor_dist

[tmp beacon_idx] = sort( beacon_vectors, 2);
C_k = beacon_idx(:,1:k);

%now calculate the routing path matrix
if p==0
  routing_path = cell(N,N);
  routing_hops = zeros(N,N); 
  %routing_dist = zeros(N,N);
else
  routing_path = cell(p,1);
  routing_hops = zeros(p,1); 
  %routing_dist = zeros(p,1);
end  
routing_load = zeros(N,1);

%now pre-compute delta(from, to,i) as \delta_i(from,to) in the
%paper
%fprintf('BVR: computing delta function.\n');
if 0 %substitute by the C implementation
delta = zeros(N,N, k);
for dest = 1:N
  C_k_dest = C_k(dest,:);
  for src = 1:N
    delta_plus = 0;
    delta_minus = 0;
    for i = 1:k
      idx = C_k_dest(i);
      diff = beacon_vectors(src,idx) - beacon_vectors(dest, idx);
      if diff>0
	delta_plus = delta_plus + diff;
      else
	delta_minus = delta_minus - diff;
      end
      delta(src, dest, i) = 10*delta_plus + delta_minus;
    end
  end
end
end

%delta = bvr_delta(beacon_vectors, C_k);
%delta = bvr_delta_new(beacon_vectors, beacon_vectors, C_k);

neighbor_2hops = cell(N,1);
neighbor_via = cell(N,1);
if fetch_2hop_info % pre-fetch it to improve speed
  fprintf('prefetching 2hop neighbor info\n');
  for i=1:N
    %neighbor_2hops{i} = [];
    %for nb = neighbor_cell{i}
    %  neighbor_2hops{i} = union(neighbor_2hops{i}, neighbor_cell{nb});
    %end
    %neighbor_2hops{i} = setdiff(neighbor_2hops{i}, [neighbor_cell{i} ...
    %		    i]);
    [dist via] = dijkstra(neighbor_hops, i);
    neighbor_2hops{i} = find(dist==2);
    tmp = find(dist~=2);
    via(tmp) = zeros(length(tmp),1);
    neighbor_via{i} = sparse(via);
  end
  fprintf('prefetching 2hop neighbor before failure update.Done\n');
end


%%%%use DHT to store node coordinates, see the overhead
%storage overhead per beacon, + packets sent
ordinary_nodes = setdiff(1:N, beacons);
storage_load = zeros(length(beacons),1);
dht_packets = 0;

if study_dht
  for i=ordinary_nodes
    b_k = consistent_hash(i,beacons,N); %b_k \in 1.. length(beacons)
    storage_load(b_k) = storage_load(b_k)+1;
    %dht_packets = dht_packets + beacon_vectors(i,b_k);
    [ tmp_path, tmp_hops, tmp_dist, tmp_load, used_2hop_info,tran_dist] ...
	= bvr_routing(beacons(b_k), i, k, debug, beacons,...
		      beacon_vectors, node_parents, C_k, ...
		      neighbor_hops, neighbor_cell, neighbor_dist, ...
		      USE_STANDARD_BVR, zeros(N,1), zeros(N,1),...
		      run_scoped_flooding, ...
		      fetch_2hop_info, neighbor_2hops, neighbor_via);
    dht_packets = dht_packets + beacon_vectors(i,b_k) + tran_dist;
  end
  if length(beacons)<10
    storage_load
  end
  
  fprintf('BVR: dht_packets:%d; mean storage load:%.2f, max:%.2f\n', ...
	  dht_packets, mean(storage_load), max(storage_load));
end

%%%%
start_t = cputime;
cur_t = cputime;
fprintf('BVR: computing routes\n');

%neighbor_cell reflects failures
for i=1:N
  neighbor_cell{i}=find(neighbor_hops(i,:) & ~failure_nodes);
end
tmp = find(failure_nodes);
neighbor_hops(tmp,:) = zeros(length(tmp), N);
neighbor_hops(:, tmp) = zeros(N, length(tmp));

if fetch_2hop_info & sum(failure_nodes)>0 % pre-fetch again due to failure
  neighbor_2hops = cell(N,1);
  neighbor_via = cell(N,1);
  fprintf('prefetching 2hop neighbor info\n');
  for i=1:N
    %neighbor_2hops{i} = [];
    %for nb = neighbor_cell{i}
    %  neighbor_2hops{i} = union(neighbor_2hops{i}, neighbor_cell{nb});
    %end
    %neighbor_2hops{i} = setdiff(neighbor_2hops{i}, [neighbor_cell{i} ...
    %		    i]);
    [dist via] = dijkstra(neighbor_hops, i);
    neighbor_2hops{i} = find(dist==2);
    tmp = find(dist~=2);
    via(tmp) = zeros(length(tmp),1);
    neighbor_via{i} = sparse(via);
  end
  fprintf('prefetching 2hop neighbor after failure update.Done\n');
end

if p == 0
  srcloop = 1:N;
  total_count = N;
else
  srcloop = test_pairs(:,1)';
  total_count = p;
end

count = 0; 
n_failed = 0;
for src=srcloop
  new_t = cputime;
  if new_t-cur_t>60
    fprintf('BVR: progress:%.1f%%. ETA %.1f sec\n',count*100/total_count, (new_t-start_t)/count*(total_count-count));
    cur_t = new_t;
  end
  if p==0
    destloop = 1:N;
  else
    destloop = [test_pairs(count+1,2)];
  end
  for dest=destloop
    count = count +1;
    [ routing_path{count}, routing_hops(count), routing_dist(count), ...
      routing_load, used_2hop_info, transmission_dist(count)] ...
	= bvr_routing(src, dest, k, debug, beacons,...
		      beacon_vectors, node_parents, C_k, ...
		      neighbor_hops, neighbor_cell, neighbor_dist, ...
		      USE_STANDARD_BVR, failure_nodes, routing_load, ...
		      run_scoped_flooding, fetch_2hop_info, ...
		      neighbor_2hops, neighbor_via);
    if fetch_2hop_info
      nodes_with_2hop_info(used_2hop_info) = ones(length(used_2hop_info), ...
					    1);
    end
    %routing_load(routing_path{count}) =
    %routing_load(routing_path{count})+1;
    %now routing_load is updated within the bvr_routing
    if routing_hops(count)==Inf
      n_failed = n_failed +1;
    end

  end
end

failure_rate = n_failed/count;
fprintf('BVR: failure rate:%.2f%%\n', failure_rate*100);

for i=1:N
  state_per_node(i,1) = sum(neighbor_dist(i,:))+length(neighbor_2hops{i})*nodes_with_2hop_info(i);
end

fprintf('%.3f node used 2 hop information\n', sum(nodes_with_2hop_info)/N);

if fetch_2hop_info
  for i=1:N
    if nodes_with_2hop_info(i)
      neighbor_i = find(neighbor_dist(i,:));
      tmp = neighbor_dist(neighbor_i,:);
      entries = sum(tmp(:));
      extra_control_traffic_2h(i) = entries*length(beacons);
    end
  end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function [routing_path, routing_hops, routing_dist, routing_load, used_2hop_info, transmission_dist] = ...
    bvr_routing(src,dest,k,debug, beacons,...
    beacon_vectors, node_parents, C_k, ...
    neighbor_hops, neighbor_cell, neighbor_dist, USE_STANDARD_BVR, ...
		failure_nodes, routing_load, run_scoped_flooding, ...
		fetch_2hop_info, neighbor_2hops, neighbor_via)
%init packet header
%packet_dst = beacon_vectors(dest,C_k(dest,:));
routing_dist = 0;
transmission_dist = 0;
used_2hop_info = [];
used_count = 0;
if failure_nodes(src) | failure_nodes(dest)
  %fprintf('source or dest failure\n');
  transmission_dist = Inf;
  routing_dist = Inf;
  routing_path = [];
  routing_hops = Inf;
  return
end
packet_delta_min = repmat(inf,k,1);
curr = src;
routing_path = [src]; routing_hops = 0; 

%routing_load(src) = routing_load(src)+1;
%start routing
if debug
  fprintf('%d->%d:',src,dest);
end

while curr~=dest
  %now the packet arrives at node curr
  if neighbor_hops(curr,dest)>0 %not described in BVR!
    next = dest;
  else
    %first update packet header
    delta_curr =bvr_delta_new(beacon_vectors(curr,:), ...
			      beacon_vectors(dest,:), C_k(dest, ...
						  :),USE_STANDARD_BVR); 
	
    packet_delta_min = min(packet_delta_min, delta_curr(:));
    %try greedy forwarding first
    next = greedy_forward(k, beacon_vectors, dest, C_k, USE_STANDARD_BVR, ...
			  neighbor_cell{curr}, packet_delta_min);
    
    if (next<0) & fetch_2hop_info %on demand 2 hop info fetch
      if debug 
	fprintf('on demand fetching 2 hop neighbors\n');
      end
      used_count = used_count +1;
      used_2hop_info(used_count) = curr;
      next = greedy_forward(k, beacon_vectors, dest, C_k, USE_STANDARD_BVR, ...
			    neighbor_2hops{curr}, ...
			    packet_delta_min);
      if next>0 %it's second hop next, not really the real next
	next = neighbor_via{curr}(next);
      end
    end
    
    if next<0 %greedy failed, use fallback mode
      fallback_bcn = C_k(dest,1);
      if beacon_vectors(curr, fallback_bcn)==Inf
	%not connected
	routing_hops = Inf;
	transmission_dist = Inf;
	routing_dist = Inf;
	fprintf('not connected to the fall back beacon\n');
	break
      elseif beacons(fallback_bcn)~=curr 
	if debug
	  fprintf(' FB(%d)',beacons(fallback_bcn));
	end
	next = node_parents(curr, fallback_bcn);
	%failure handling
	if failure_nodes(next)>0
	  routing_hops = Inf;
	  transmission_dist = Inf;
	  routing_dist = Inf;
	  %fprintf(' parent to beacon is dead\n');
	  if debug
	    fprintf(' dead\n');
	  end
	  break
	end
      else
	%totally failed, have to broadcast. 
	if debug
	  fprintf(' FAILED ');
	  fprintf('broadcasting within %d\n',beacon_vectors(dest, fallback_bcn));
	end
	routing_hops = Inf;
	%scoped flooding
	if run_scoped_flooding
	  %[b_snd, b_rcv, nexthop, label] = BFS(curr, beacon_vectors(dest, fallback_bcn), ...
	  %				       neighbor_dist,1,-1);
	  %%%%%%%%don't use above--deprecated
	  [found, bd_nodes] = scoped_flooding( curr, beacon_vectors(dest, ...
						  fallback_bcn), ...
					       neighbor_hops, dest);
	  
%	  routing_path = [ routing_path find(nexthop) dest];
% $$$ 	  if nexthop(dest)<=0
% $$$ 	    nexthop
% $$$ 	    routing_hops = Inf;
% $$$ 	    routing_dist = Inf;
% $$$ 	    fprintf('BFS failed\n');
%	  else
%	    routing_dist = routing_dist + b_snd;
%	  end

	  routing_load(bd_nodes) = routing_load(bd_nodes)+1;
	  transmission_dist = transmission_dist + length(bd_nodes);
	  routing_dist = routing_dist + ...
	      beacon_vectors(dest, fallback_bcn);
	  if ~found
	    routing_hops=Inf;
	    transmission_dist = Inf;
	    routing_dist = Inf;
	    fprintf('scoped flooding failed\n');
	  else
	    routing_path = [routing_path bd_nodes dest];
	    %routing_load(dest) = routing_load(dest) + 1;
	  end
	else
	  %don't simulate scoped flooding. indicate a failure only
	  routing_hops = Inf;
	  transmission_dist = Inf;
	  routing_dist = Inf;
	end
	break
      end
    end %end of if greedy failed
  end
  if debug
    fprintf(' %d ',next);
  end
  %      if routing_hops(count)>50
  %	fprintf('src:%d,dest:%d:%d,%d\n',src,dest,routing_hops(count),next);
  %      end
  routing_hops = routing_hops +1;
  routing_path(routing_hops) = next;
  routing_load(curr) = routing_load(curr)+1;
  transmission_dist = transmission_dist + neighbor_dist(curr, next);
  routing_dist = routing_dist + neighbor_dist(curr, next);
  curr = next;
  
end
if debug
  fprintf('\n');
end

function next = greedy_forward(k, beacon_vectors, dest, C_k, USE_STANDARD_BVR, ...
			  neighbors, packet_delta_min);
next = -1;
neighbor_delta = ...
    bvr_delta_new(beacon_vectors(neighbors,:), ...
		  beacon_vectors(dest,:), C_k(dest,:), USE_STANDARD_BVR);
%	sqrt(neighbor_delta) == ...
%	    euclidean_distance(beacon_vectors(neighbor_cell{curr},C_k(dest,:)), ...
%			       beacon_vectors(dest,C_k(dest,:)))

if USE_STANDARD_BVR
  kloop = k:-1:1;
else
  kloop = 1:1;
end
if length(neighbors)>0
  for i=kloop
    [v next_candidate] = min( neighbor_delta(:,:,i) ); %%
    if v<packet_delta_min(i)
      %next = next_candidate;
      next = neighbors(next_candidate);
      break
    end
  end
end