function [K_nn]= reducematrix( K_gg, G_mn, pv_m_set, pv_n_set, msize, flag );
% REDUCEMATRIX [K_nn]= reducematrix( K_gg, G_mn, pv_m_set, pv_n_set, msize, flag )
% Reduces K_gg to K_nn
%          
%          K_nn^~   K_nm
% K_gg = [               ]
%          K_mn     K_mm
%
% if flag == 1
%  	K_nn = K_nn^~ + K_nm * G_mn
% else
% K_nn = K_nn^~ + K_nm * G_mn + G_mn^T * K_mn + G_mn^T * K_mm * G_mn
%
% m dof to be condensed 
% n dof to be retained
% g total dof g = m + n
%
% Input:  K_gg    	matrix (g x g)  
%		  G_mn      Reduction matrix from constraints (m x n)
%         pv_m_set  partitioning vector for the m dof wrt g 
%         pv_n_set  partitioning vector for the n dof wrt g
%		  msize		number of m-set dof
%         flag      flag for type of condensation, see above 
%
% Output: K_nn      reduced matrix (n x n) 
%
% Called by  Reducematrixfromgtof (mf), Normalmodes (ms)
% Calls none 
% where ms= m-script, mf= m-function, mex= executable function (c-code)
%
% pv_m_set and pv_n_set index vectors are with respect to the g- set.
%
% All input and output variables are assumed to be in the workspace.

if msize > 0 
		
%   Partition K_gg

	K_mm = K_gg( pv_m_set, pv_m_set);
	K_mn = K_gg( pv_m_set, pv_n_set);
	K_nm = K_gg( pv_n_set, pv_m_set);
	K_nn = K_gg( pv_n_set, pv_n_set);
		
%   Reduce K_gg to K_nn

	K_nn = K_nn + K_nm * G_mn;
	
	if  flag ~= 1
	
		K_mn = K_mn + K_mm * G_mn;
		K_nn = K_nn + G_mn' * K_mn;
	
	end
	
else     
	
% 	msize=0 just return input

	K_nn = K_gg;
	
end
