 
 			// Loop through each pair of labels
 			T minCost = grid->getLabelingCost(result);
-			for (int alpha = 1; alpha < grid->getNumLabels(); alpha++) {
-				for (int beta = 0; beta < alpha; beta++) {
+			for (int alpha = 0; alpha < grid->getNumLabels(); alpha++) {
+				printf("%lf,%d\n", TIMER_END, grid->getLabelingCost(result));
+			
+				for (int beta = alpha + 1; beta < grid->getNumLabels(); beta++) {
 					// Now, construct the graph with alpha and beta as terminal nodes
 					gpuGraphCutConstructGraphKernel<<<gridSize, blockSize>>>(dGrid, dResult, dHeight, dExcess, dResidue, alpha, beta);
 
 					// We need to keep this data synced between CPU and GPU
 					// every time we perform global relabeling
 					gpu.copyRawMat(&height, dHeight);
-					gpu.copyRawMat(&excess, dExcess);
 					gpu.copyRawMat(&residue, dResidue);
 
 					// Run push-relabel until convergence
 					int iters = 0;
 					while (true) {
+						// Mark each tile as active or inactive
+						gpu.clearMat(dActive, isActive.getTotalSize());
+						gpuFindActiveTiles<<<gridSize, blockSize>>>(dGrid, dHeight, dExcess, dActive);
+						gpu.copyRawMat(&isActive, dActive);
+						
+						// Build a list of active tiles
+						if (buildActiveTileList(&isActive, activeList) == 0) break;
+						
 						// Global relabeling optimization on the CPU
-						if (graphCutGlobalRelabel(grid, &height, &excess, &residue)) {
-							// No more active nodes, we're done
-							break;
-						}
+						graphCutGlobalRelabel(grid, &height, &residue);
 						
 						// Copy relabeled height data back to the GPU
 						gpu.copyRawMat(dHeight, &height);
 
 						// Run normal push and relabel on the GPU
 						for (int i = 0; i < GINF_GRAPHCUT_NUM_ITERS_PER_GLOBAL_RELABEL; i++) {
-							gpuGraphCutPushKernel<<<gridSize, blockSize>>>(iters & 1, dGrid, dHeight, dExcess, dResidue);
-							cudaThreadSynchronize();
+							gpuGraphCutPushKernel<<<gridSize, blockSize>>>(dGrid, dHeight, dExcess, dResidue);
+							gpuGraphCutRelabelKernel<<<gridSize, blockSize>>>(dGrid, dHeight, dExcess, dResidue);
 							
-							gpuGraphCutRelabelKernel<<<gridSize, blockSize>>>(iters & 1, dGrid, dHeight, dExcess, dResidue);
-							cudaThreadSynchronize();
 							iters++;
-							
 						}
 						
 						// Copy graph data back to CPU for global relabeling
 						gpu.copyRawMat(&height, dHeight);
-						gpu.copyRawMat(&excess, dExcess);
 						gpu.copyRawMat(&residue, dResidue);
 					}
 

	}
	
	template <typename T>
	__global__ void gpuGraphCutPushKernel(bool iter, GpuGrid<T> *grid, int *height, T *excess, T *residue) {
		GINF_DECL_ALT_X(iter); GINF_DECL_Y;
		MatDim dim(grid->getWidth(), grid->getHeight());

		if (dim.isValid(x, y) && excess[dim.idx(x, y)] && height[dim.idx(x, y)] < grid->getNumNodes() + 2) {

		// Start with the intial labeling
		result->copyFrom(initial);
		
		// Prepare block sizes
		dim3 blockSize(GINF_GRAPHCUT_BLOCK_WIDTH, GINF_GRAPHCUT_BLOCK_HEIGHT);
		dim3 gridSize(GINF_DIM(grid->getWidth(), blockSize.x), GINF_DIM(grid->getHeight(), blockSize.y));
		
		// Boolean matrix indicating whether a node is reachable from the source in the residue graph
		Matrix<int> reachable(2, grid->getWidth(), grid->getHeight());
		
		// Boolean matrix indicate whether a block contains an active node or not
		Matrix<int> isActive(2, gridSize.x, gridSize.y);
		
		// List of active blocks (xy coordinates)
		int2 activeList[gridSize.x * gridSize.y];

		// Create a temp matrix to store new labelings
		Matrix<int> newResult(2, grid->getWidth(), grid->getHeight());
