Index: /issm/trunk/src/m/solvers/solver_adjoint_linear.m
===================================================================
--- /issm/trunk/src/m/solvers/solver_adjoint_linear.m	(revision 5998)
+++ /issm/trunk/src/m/solvers/solver_adjoint_linear.m	(revision 5999)
@@ -5,13 +5,18 @@
 %      femmodel =solver_adjoint_linear(femmodel)
 
-	[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
-	
-	[K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
-	%K_ff=transpose(K_ff);
-	p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters,true);
+	%Get parameters
+	kffpartitioning=femmodel.parameters.Kff;
+
+	if kffpartitioning,
+		[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
+		p_f = Reduceload( p_f, K_fs, femmodel.ys);
+	else
+		[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
+		[K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
+		p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters,true);
+	end
 
 	displaystring(femmodel.parameters.Verbose>1,'%s%g','      condition number of stiffness matrix: ',condest(K_ff));
 	u_f=Solver(K_ff,p_f,[],femmodel.parameters);
-	
 	u_g= Mergesolutionfromftog( u_f, femmodel.ys, femmodel.nodesets,femmodel.parameters,true); 
 
Index: /issm/trunk/src/m/solvers/solver_diagnostic_nonlinear.m
===================================================================
--- /issm/trunk/src/m/solvers/solver_diagnostic_nonlinear.m	(revision 5998)
+++ /issm/trunk/src/m/solvers/solver_diagnostic_nonlinear.m	(revision 5999)
@@ -29,9 +29,7 @@
 
 		if kffpartitioning,
-			disp('-------------- file: solver_diagnostic_nonlinear.m line: 31'); 
 			[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters);
 			p_f = Reduceload( p_f, K_fs, femmodel.ys);
 		else
-			disp('-------------- file: solver_diagnostic_nonlinear.m line: 35'); 
 			[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters);
 			[K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
Index: /issm/trunk/src/m/solvers/solver_linear.m
===================================================================
--- /issm/trunk/src/m/solvers/solver_linear.m	(revision 5998)
+++ /issm/trunk/src/m/solvers/solver_linear.m	(revision 5999)
@@ -5,12 +5,18 @@
 %      femmodel =solver_linear(femmodel)
 
-	[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
-	
-	[K_ff, K_fs] = Reducematrixfromgtof( K_gg,  femmodel.nodesets,femmodel.parameters); 
-	p_f = Reduceloadfromgtof( p_g,  K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters);
+	%Get parameters
+	kffpartitioning=femmodel.parameters.Kff;
+
+	if kffpartitioning,
+		[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
+		p_f = Reduceload( p_f, K_fs, femmodel.ys);
+	else
+		[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
+		[K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
+		p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters);
+	end
 
 	displaystring(femmodel.parameters.Verbose>1,'%s%g','      condition number of stiffness matrix: ',condest(K_ff));
 	u_f=Solver(K_ff,p_f,[],femmodel.parameters);
-	
 	u_g= Mergesolutionfromftog( u_f, femmodel.ys, femmodel.nodesets,femmodel.parameters); 
 
Index: /issm/trunk/src/m/solvers/solver_stokescoupling_nonlinear.m
===================================================================
--- /issm/trunk/src/m/solvers/solver_stokescoupling_nonlinear.m	(revision 5998)
+++ /issm/trunk/src/m/solvers/solver_stokescoupling_nonlinear.m	(revision 5999)
@@ -4,4 +4,7 @@
 %   Usage:
 %      [femmodel]=solver_couplingstokes_nonlinear(femmodel,conserve_loads)
+
+	%Get parameters
+	kffpartitioning=femmodel.parameters.Kff;
 
 	%initialize solution vector
@@ -31,11 +34,14 @@
 		old_uf_horiz=uf_horiz;
 
-		[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
-
-		[K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
-		p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters);
+		if kffpartitioning,
+			[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
+			p_f = Reduceload( p_f, K_fs, femmodel.ys);
+		else
+			[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
+			[K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
+			p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters);
+		end
 
 		uf_horiz=Solver(K_ff,p_f,old_uf_horiz,femmodel.parameters);
-
 		ug_horiz= Mergesolutionfromftog( uf_horiz, femmodel.ys, femmodel.nodesets,femmodel.parameters); 
 
@@ -48,8 +54,12 @@
 %		%Then compute vertical velocity
 %		femmodel=SetCurrentConfiguration(femmodel,DiagnosticVertAnalysisEnum);
+%	if kffpartitioning,
 %		[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
-%
-%		[K_ff, K_fs] = Reducematrixfromgtof( K_gg,  femmodel.nodesets,femmodel.parameters); 
-%		p_f = Reduceloadfromgtof( p_g,  K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters);
+%		p_f = Reduceload( p_f, K_fs, femmodel.ys);
+%	else
+%		[K_gg,K_ff,K_fs,p_g,p_f,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
+%		[K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
+%		p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters);
+%	end
 %
 %		displaystring(femmodel.parameters.Verbose>1,'%s%g','      condition number of stiffness matrix: ',condest(K_ff));
@@ -61,4 +71,3 @@
 
 	end
-
 end
Index: /issm/trunk/src/m/solvers/solver_thermal_nonlinear.m
===================================================================
--- /issm/trunk/src/m/solvers/solver_thermal_nonlinear.m	(revision 5998)
+++ /issm/trunk/src/m/solvers/solver_thermal_nonlinear.m	(revision 5999)
@@ -8,4 +8,7 @@
 	count=1;
 	converged=0;
+
+	%Get parameters
+	kffpartitioning=femmodel.parameters.Kff;
 
 	%stiffness and load generation only:
@@ -20,8 +23,12 @@
 	while(~converged),
 
-		[K_gg,K_ff,K_fs,p_g,p_f,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
-
-		[K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
-		p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters);
+		if kffpartitioning,
+			[K_gg,K_ff,K_fs,p_g,p_f,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
+			p_f = Reduceload( p_f, K_fs, femmodel.ys);
+		else
+			[K_gg,K_ff,K_fs,p_g,p_f,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
+			[K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets,femmodel.parameters); 
+			p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,femmodel.parameters);
+		end
 
 		displaystring(femmodel.parameters.Verbose,'%s%g','   condition number of stiffness matrix: ',condest(K_ff));