function [A1, B1, C1, const] = fixrs1( A, B, C, r, s);
  
% QAP instance with (A,B,C) is given. A,B symmetric assumed
% new instance is generated in which r is assigned to s
% resulting in a problem of smaller dimension + constant
% one 'step' down!
% call: [A1, B1, C1, const] = fixrs1( A, B, C, r, s);

  
n = size(A,1);
mi = min( r, s);
ma = max( r, s);
if (mi<1) | (ma > n);
disp(' input wrong: r or s outside limits. nothing done.');
A1=[]; B1=A1; C1=A1; const=A1;
return
end;

N = 1:n;
I = setdiff( N, r);   % take away element r
J = setdiff( N, s);   % take away element s
const = A(r,r) * B(s,s) + C(r,s);   % contribution from assigning r to s
A1 = A(I,I);          % delete row and column r 
B1 = B(J,J);          % delete row and column s
ar = A(:, r);         % extract column r 
ar = ar( I);          % ignore main diag element A(r,r)
bs = B(:, s);         % extract column s
bs = bs( J);          % ignore main diag element B(s,s)
C1 = C(I,J) + 2 * ar*bs';  % remove row r and col s from C
                      % and add 2 the contribution from partial ass.
		      % here symmetry of A and B is assumed
		      

		      
% for testing 
% I=[2 5 4 7 9 1 3 6 8 12 11 10];
% pom = eye(n); X=pom(I,:);
% test1= trace( A*X*B*X' + C*X');
% N = 1:n; II=setdiff( N, r); JJ=setdiff( N, s);  X1 = X(II,JJ);
% [A1, B1, C1, const] = fixrs1( A, B, C, r, s);
% test2= trace( A1*X1*B1*X1' + C1*X1') + const; norm(test1 - test2)
