science-ation/admin/anneal.inc.php

<?php

class annealer {

	var $num_buckets;
	var $bucket;
	var $bucket_cost;
	var $cost;
	var $start_temp, $start_moves;
	var $cost_function_callback;
	var $iterations;
	var $items_per_bucket;
	var $rate;
	
	function annealer($num_buckets, $start_temp, $start_moves, $rate,
				$cost_function_cb, $items)
	{
		$this->num_buckets = $num_buckets;
		$this->start_temp = $start_temp;
		$this->start_moves = $start_moves;
		$this->cost_function_callback = $cost_function_cb;
	
		$this->bucket_cost = array();
		$this->bucket = array();
		$this->iterations = 0;
		$this->items_per_bucket = count($items) / $num_buckets;
		$this->rate = $rate;
		$ipb = ceil($this->items_per_bucket);
		$i=0;
		$this->cost = 0;
		for($x=0; $x<$num_buckets; $x++) {
			unset($b);
			$b = array();
			for($y=0;$y<$ipb; $y++) {
				if($i == count($items)) break;
				$b[] = $items[$i];
				$i++;
			}
			$this->bucket[] = $b;
			$c = $this->cost_function($x);
			$this->bucket_cost[] = $c;
			$this->cost += $c;
		}
		TRACE("Annealer setup: T={$this->start_temp}, ".
				"M={$this->start_moves}, Bkts={$this->num_buckets}, ".
				"Cost={$this->cost}\n");
	}


	function pick_move()
	{
		/* Pick a bucket and item */
		while(1) {
			$b1 = rand(0, $this->num_buckets - 1);
			if(count($this->bucket[$b1]) > 0) break;
		}
		$i1 = rand(0, count($this->bucket[$b1]) -1);

		/* Pick a csecond bucket that is different thatn the first */
		$b2 = rand(0, $this->num_buckets - 2);
		if($b2 >= $b1) $b2++;
	
		/* Picket an item, or a blank, in the second bucket */
		$i2 = rand(0, count($this->bucket[$b2]));
		if($i2 == count($this->bucket[$b2])) $i2 = -1;
//		TRACE("Move ($b1,$i1)<->($b2,$i2)\n");
		return array($b1, $i1, $b2, $i2);
	}
	
	function cost_function($b)
	{
		$bkt = $this->bucket[$b];
		$cb = $this->cost_function_callback;
		$c = $cb($this, $b, $bkt);
//		$this->print_bucket($b);
//		print("Computed cost to be: $c\n");
		return $c;
	}
	
	function compute_delta_cost($move) 
	{
		list($b1, $i1, $b2, $i2) = $move;

		$cost = 0;

		$b1_old = $this->bucket[$b1];
		$b2_old = $this->bucket[$b2];
		
		$b1_new = array();
		$b2_new = array();
		/* Make 2 new bucket lists */
		for($x=0; $x<count($b1_old); $x++) {
			$id = $b1_old[$x];
			if($x == $i1) {
				/* Swap or remove this index */
				if($i2 != -1) $b1_new[] = $b2_old[$i2];
			} else {
				$b1_new[] = $id;
			}
		}
		for($x=0; $x<count($b2_old); $x++) {
			$id = $b2_old[$x];
			if($x == $i2) {
				/* Swap or remove this index */
				$b2_new[] = $b1_old[$i1];
			} else {
				$b2_new[] = $id;
			}
		}
		if($i2 == -1) $b2_new[] = $b1_old[$i1];
		

		/* Assign the new item lists to the buckets */
		$this->bucket[$b1] = $b1_new;
		$this->bucket[$b2] = $b2_new;

		/* Compute costs */
		$cost -= $this->bucket_cost[$b1];
		$cost -= $this->bucket_cost[$b2];

		$c1 = $this->cost_function($b1);
		$c2 = $this->cost_function($b2);

		$cost += $c1 + $c2;

		/* Return to the original bucket lists */
		$this->bucket[$b1] = $b1_old;
		$this->bucket[$b2] = $b2_old;
	
		return array($cost, array($c1, $b1_new, $c2, $b2_new));
	}

	function accept_move($move, $movedata)
	{
		list($b1, $i1, $b2, $i2) = $move;
		list($c1, $b1_new, $c2, $b2_new) = $movedata;
//		TRACE("Old buckets:\n");
//		$this->print_bucket($b1);
//		$this->print_bucket($b2);
		$this->bucket[$b1] = $b1_new;
		$this->bucket[$b2] = $b2_new;
//		TRACE("New buckets:\n");
//		$this->print_bucket($b1);
//		$this->print_bucket($b2);
		$this->bucket_cost[$b1] = $c1;
		$this->bucket_cost[$b2] = $c2;
	}

	function anneal()
	{
		$temperature = $this->start_temp;
		$current_cost = $this->cost;
		$last_cost = 0;
		$last_cost_count = 0;

		if($this->num_buckets <= 1) {
			TRACE("Only one Bucket, nothing to anneal.\n");
			return;
		}
//		$this->print_buckets();
		while(1) {
			$moves = $this->start_moves;
			for($m = 0; $m<$moves; $m++) {
//				$this->print_buckets();
				/* Pick 2 moves at random */
				$move = $this->pick_move();
				/* See what the new cost is compared to the old */
				list($delta_c, $movedata) = 
					$this->compute_delta_cost($move);


				$r = floatval(rand()) / floatval(getrandmax());
				/* Decide if we want to keep it */
				$e = exp(-$delta_c / $temperature);
//				TRACE("r=$r, exp=$e, delta=$delta_c\n");
				if($r < exp(-$delta_c / $temperature)) {
					/* Yes, we do, record the move */
					$this->accept_move($move, $movedata);
					$current_cost += $delta_c;
					$n_accepted++;
					if($current_cost < $this->cost)
						$this->cost = $current_cost;
//					TRACE("Move accepted, cost=$current_cost\n");
				} else {
//					TRACE("Move rejected\n");
				}
				$this->iterations++;
				if($this->iterations % 50000 == 0) {
					TRACE("   {$this->iterations} iterations, cost={$this->cost}, temperature=$temperature\n");
//					$this->print_buckets();
				}
			}
			if($this->cost == $last_cost) {
				$last_cost_count ++;
			} else {
				$last_cost = $this->cost;
				$last_cost_count=0;
			}
			
			if($temperature < 0.1 && $last_cost_count > 10) 
				break;
//			TRACE("Cost is {$this->cost}\n");
			$temperature *= $this->rate;
		}
		TRACE("Annealing complete.  {$this->iterations} iterations.  Final cost is {$this->cost}\n");
	}

	function print_bucket($x)
	{
		$b = $this->bucket[$x];
		print("Bucket $x: (cost: {$this->bucket_cost[$x]})\n");
		print("   ");
		for($y=0;$y<count($b); $y++) {
			print("{$b[$y]} ");
		}
		print("\n");
	}
	function print_buckets()
	{
		for($x=0; $x<$this->num_buckets; $x++) {
			$this->print_bucket($x);
		}
	}
}
?>