% hopf.m for showing the biproduct pairs

  fprintf('hopf.m --- Match the single indices (left) %s\n',... 
            'with double index pairs (right)')
for n=[4 7] 
  m=(n-1)*n/2; Index_Pair=[]; 
  L='--------------------------------------------------------------';L=[L L];
  S='                                                              ';S=[S S];

for I=1:m
   i = (-1+sqrt(1+8*I))/2;  % Positive root of (i+1)*i/2 = I
  ii = floor(i); i2 = ceil(i);  ni = (ii+1)*ii/2;  % 1st diag
  idiff = i2 - ii; 
  if idiff< 1 
     is = ii+1;  js = ii;   % diag case
  else
     is = ii+2;  js = I-ni; % from 1st diag
  end
  Index_Pair = [Index_Pair  [is; js]];
end
  for i=2:n
    i1=(i-2)*(i-1)/2+1; i2=(i-1)*i/2; 
    fprintf('%3d ',i1:i2); fprintf('%s',S(1:4*(n-i))); 
    fprintf('(%2d,%2d) ',Index_Pair(:,i1:i2));  fprintf('\n')
  end
    fprintf('%s\n',L(1:4*(n-1)))
end
    fprintf('\t\t%%=================Reverse the biproduct index formulae...\n')
  fprintf('hopf.m --- Match the double index pairs (left) %s\n',...
                  'with single indices (right)') 
for n=[4 7] 
 for p=2:n
     U=[];
  for q=1:p-1
   u = (p-2)*(p-1)/2 + q; U=[U u]; fprintf('(%2d,%2d) ',p,q);
  end
    fprintf('%s',S(1:8*(n-p))); fprintf('%3d ',U); fprintf('\n')
 end
    fprintf('%s\n',L(1:8*(n-1)))
end
