  function [fitted,vaf,weights,end_condition] = partitionfit_twomode(proxtm,member)

% PARTITIONFIT_TWOMODE provides a least-squares approximation to a two-mode
% proximity matrix based on a given collection of partitions.
%
% syntax: [fitted,vaf,weights,end_condition] = partitionfit_twomode(prox,member)
%
% PROXTM is the nrow x ncol two-mode input proximity matrix (with a zero main diagonal
% and a dissimilarity interpretation); MEMBER is the m x n (n = nrow+ ncol) matrix
% indicating cluster membership, where each row corresponds to a specific
% partition (there are m partitions in general); the columns of MEMBER
% are in the same order used for PROXTM, with rows first and columns to follow.
% FITTED is an nrow x ncol matrix fitted to PROXTM (through least-squares)
% constructed from the nonnegative weights given in the m x 1 WEIGHTS
% vector corresponding to each of the partitions.  VAF is the variance-
% accounted-for in the proximity matrix PROXTM by the fitted matrix FITTED.
% END_CONDITION should be zero for a normal termination of the optimization
% process.

[nrow,ncol] = size(proxtm);
n = nrow + ncol;
nsize = nrow*ncol;
npart = size(member,1);



partvec = zeros(npart,nsize);
datavec = zeros(1,nsize);
fitted = zeros(nrow,ncol);
proxave = zeros(nrow,ncol);

index = 0;

for i = 1:nrow
    for j = 1:ncol
        
      
        index = index + 1;        
        datavec(index) = proxtm(i,j);
        
        for k = 1:npart
            
            if(member(k,i) ~= member(k,(j+nrow)))
                 partvec(k,index) = 1;
            end
            
        end
    end
end
    
identity = eye(nsize);
equality_flag = zeros(npart,1);
constraint_constant = zeros(npart,1);

[solution,kuhn_tucker,iterations,end_condition] = ...
    dykstra(datavec,identity,partvec,constraint_constant,equality_flag);

weights = (.5)*kuhn_tucker;

fittedvec = datavec - solution';

index = 0;

for i = 1:nrow
    for j = 1:ncol
        
        index = index + 1;
        fitted(i,j) = fittedvec(index);
    
        
    end
end

aveprox = sum(sum(proxtm))/nsize;

for i = 1:nrow
    for j = 1:ncol
        
         proxave(i,j) = aveprox;
        
    end
end

diff = sum(sum((proxtm - fitted).^2));
denom = sum(sum((proxtm - proxave).^2));

vaf = 1 - (diff/denom);