function W = bnd_share_matrix(bnd,eps)

%BND_SHARE_MATRIX computes a boundary share matrix from a boundary in ".bnd" format 
%with row-number identifier. The output is row-normalized form of this matrix. 

% NOTES: (1) If POLYFORM gives an error message, try BND_SHARE_MATRIX_QUIET
%     
%        (2) If your final matrix shows zero rows for polygons that clearly share f
%            faces with other polygons, then the default "eps" parameter in the subroutine, 
%            BND_SHARES, may be too small. Try a larger value.


%Written by: TONY E. SMITH, 1/4/02
%Modified: 3/2/04

%INPUT: bnd = (n:2) matrix describing boundaries of k polygons in the form:
%
%          [1   n1
%           x11 y11 
%           x12 y12
%           .
%           .
%           x11 y11
%           2   n2
%           x21 x22
%           .
%           .
%           .
%           k   nk
%           xk1 yk1
%           .
%           .
%           xk1 yk1]
%
%       eps = (optional) value for the error parameter in BND_SHARES

%OUTPUT: W = a full-form row-normalized matrix of boundary shares

%PROGRAMS CALLED: bnd_shares.m, row_norm.m


%Compute Shared Boundary Lengths

Z = polyform(bnd);

if nargin == 2
    
    [I,J,L] = bnd_shares(bnd,Z,eps);
    
else
    
    [I,J,L] = bnd_shares(bnd,Z);
    
end   


%Make Sparse Matrix

MAT = [I,J,L];

%Convert to Full Matrix

W = spconvert(MAT);

%************************************************************
%NOTE: Here we need to check to be sure that max of I and J is equal
%      to the number of polygons, k. If not, then we must add some zero
%      rows and columns to W, representing the ISLANDS missed.
%**************************************************************

kk = size(W,1); %size of W

Z = polyform(bnd);

k = size(Z,1); %number of polygons

if kk < k
   
   W_temp = zeros(k,k);
   W_temp(1:kk,1:kk) = W;
   W = W_temp;
   
end
   

%Row Normalize

W = full(row_norm(W));