#!/usr/bin/php
<?php

define('MAX_POWER', 10);
define('MAX_DISTANCE', 999);
ini_set('memory_limit', '256M');

# ---

function bcpow2($x) {
	return bcmul2($x,$x);
}

function bcpow3($x) {
	return bcmul3($x,$x,$x);
}

function bcmul3($a,$b,$c) {
	return bcmul2($a,bcmul2($b,$c));
}

function bcsub4($a,$b,$c,$d) {
	return bcsub2(bcsub2(bcsub2($a,$b),$c),$d);
}

function bcabs($x) {
	return str_replace('-', '', $x);
}

function bcadd2($a,$b) {
	return gmp_add($a,$b);
}

function bcsub2($a,$b) {
	return gmp_sub($a,$b);
}

function bcmul2($a,$b) {
	return gmp_mul($a,$b);
}

function bicomplex_mul($a, $b) {

	# According to https://link.springer.com/chapter/10.1007/978-3-319-24868-4_2
	# "bicomplex number is a complex number with complex coefficients."
	// Z = z1 + jz2  where z1 and z2 are complex numbers
	// Operations: Z+W := (z1+w1) + j(z2+w2)
	//             Z*W := (z1+jz2)(w1+jw2)=(z1w1 - z2w2) + j(z1w2 + z2w1)
	// Z*W := (z1+jz2)(w1+jw2) = (z1w1 - z2w2) + j(z1w2 + z2w1)
	// Z*W := ((a+bi)+j(c+di))((e+fi)+j(g+hi)) = ((a+bi)(e+fi) - (c+di)(g+hi)) + j((a+bi)(g+hi) + (c+di)(e+fi))
	// Z*W := (a+bi+cj+dij)(e+fi+gj+hij) = ((ae-bf)+i(af+be) - ((cg-dh)+i(ch+dg)) ) + j( (ag-bh)+i(ah+bg) + (ce-df)+i(cf+de) )
	// Z*W := (a+bi+cj+dij)(e+fi+gj+hij) = ae-bf+afi+bei -(cg-dh)-i(ch+dg) + j( ag-bh+ahi+bgi + ce-df+cfi+dei )
	// Z*W := (a+bi+cj+dij)(e+fi+gj+hij) = ae-bf+afi+bei -cg+dh-chi-dgi + j( ag-bh+ahi+bgi + ce-df+cfi+dei )
	// Z*W := (a+bi+cj+dij)(e+fi+gj+hij) = ae-bf+afi+bei -cg+dh-chi-dgi +  agj-bhj+ahij+bgij + cej-dfj+cfij+deij
	// Z*W := (a+bi+cj+dij)(e+fi+gj+hij) = ae-bf-cg+dh + afi+bei-chi-dgi + agj-bhj+cej-dfj + ahij+bgij+cfij+deij
	// Z*W := (a+bi+cj+dij)(e+fi+gj+hij) = (ae-bf-cg+dh) + i(af+be-ch-dg) + j(ag-bh+ce-df) + ij(ah+bg+cf+de)
	// Applying k=ij:
	// Z*W := (a+bi+cj+dk)(e+fi+gj+hk) = (ae-bf-cg+dh) + i(af+be-ch-dg) + j(ag-bh+ce-df) + k(ah+bg+cf+de)

	// Warning!!! The multiplication table on some Wikipedia articles (English,Japanese,Chinese,...)
	// is WRONG! In spanish + french Wikipedia, it is correct and also matches
	// http://3dfractals.com/docs/Article01_bicomplex.pdf as well as
	// the result of (z1w1 - z2w2) + j(z1w2 + z2w1) described in https://link.springer.com/chapter/10.1007/978-3-319-24868-4_2

	// Based on the correct multiplication table:

	// C = A*B := (a_0+a_1i+a_2j+a_3k)(b_0+b_1i+b_2j+b_3k) =
	// (a_0b_0 - a_1b_1 - a_2b_2 + a_3b_3) +
	//i(a_0b_1 + a_1b_0 - a_2b_3 - a_3b_2) +
	//j(a_0b_2 - a_1b_3 + a_2b_0 - a_3b_1) +
	//k(a_0b_3 + a_1b_2 + a_2b_1 + a_3b_0)
	// c_0 + c_1i + c_2j + c_3k

	# $c0 = $a[0]*$b[0] - $a[1]*$b[1] - $a[2]*$b[2] + $a[3]*$b[3];
	# $c1 = $a[0]*$b[1] + $a[1]*$b[0] - $a[2]*$b[3] - $a[3]*$b[2];
	# $c2 = $a[0]*$b[2] - $a[1]*$b[3] + $a[2]*$b[0] - $a[3]*$b[1];
	# $c3 = $a[0]*$b[3] + $a[1]*$b[2] + $a[2]*$b[1] + $a[3]*$b[0];
	$c0 = bcadd2(bcsub2(bcsub2(bcmul2($a[0],$b[0]),bcmul2($a[1],$b[1])),bcmul2($a[2],$b[2])),bcmul2($a[3],$b[3]));
	$c1 = bcsub2(bcsub2(bcadd2(bcmul2($a[0],$b[1]),bcmul2($a[1],$b[0])),bcmul2($a[2],$b[3])),bcmul2($a[3],$b[2]));
	$c2 = bcsub2(bcadd2(bcsub2(bcmul2($a[0],$b[2]),bcmul2($a[1],$b[3])),bcmul2($a[2],$b[0])),bcmul2($a[3],$b[1]));
	$c3 = bcadd2(bcadd2(bcadd2(bcmul2($a[0],$b[3]),bcmul2($a[1],$b[2])),bcmul2($a[2],$b[1])),bcmul2($a[3],$b[0]));

	return [ $c0, $c1, $c2, $c3 ];
}

function bicomplex_int_power($bicomplex, $n) {
	assert($n >= 1);
	$res = $bicomplex;
	for ($p=2; $p<=$n; $p++) $res = bicomplex_mul($res, $bicomplex);
	return $res;
}

function bicomplex_is_immortal($bicomplex, $power) {
	$bicomplexpow = bicomplex_int_power($bicomplex, $power);
	for ($i=0; $i<4; $i++) {
		$x1 = bcabs($bicomplex[$i]);
		$x2 = bcabs($bicomplexpow[$i]);
		if (substr($x2, -strlen($x1)) != $x1) return false;
	}
	return $bicomplexpow;
}

function _go($a, $b, $c, $d, $power, $distance, &$progr, &$found) {
	$bicomplex = [$a, $b, $c, $d];
	$bicomplexpower = bicomplex_is_immortal($bicomplex, $power);
	if ($bicomplexpower !== false) {
		$a_ = $bicomplexpower[0];
		$b_ = $bicomplexpower[1];
		$c_ = $bicomplexpower[2];
		$d_ = $bicomplexpower[3];
		// Z=z1+jz1   z1=a+bi    z2=cj+dij=cj+dk  (k=ij)
		if (($c!=0) || ($d!=0)) {
			// BICP	z1 + jz2	with z2=...+...i and z1=...+...i
			$str = "BICP ";
		} else if ($b!=0) {
			// CPLX	z1 + jz2	with z2=0+0i and z1=...+...i
			$str = "CPLX ";
		} else {
			// REAL	z1 + jz2	with z2=0+0i and z1=...+0i
			$str = "REAL ";
		}
		$str .= "| Distance=$distance | ";
		$str .= "[$a + $b i + $c j + $d k] ^ $power = [$a_ + $b_ i + $c_ j + $d_ k]";
#		if (!in_array($str,$found)) {
#			echo "$str\n";
#			$found[] = $str;
file_put_contents('/tmp/daniel_immortal_bicomplex', "$str\n", FILE_APPEND);
$found++;
#		}
	}
	#if ($progr++ % 10000 == 0) {
	#	echo "(Dist=$distance) Processing ($a, $b, $c, $d) with power $power                 \r";
	#}
}

if (file_exists(__DIR__.'/savestate.txt')) {
	$start_distance = trim(file_get_contents(__DIR__.'/savestate.txt'));
} else {
	$start_distance = 0;
}
for ($distance=$start_distance; $distance<=MAX_DISTANCE; $distance++) {
	$str = "--------- Distance: $distance ---------";
	echo "                                                                              \r";
	echo "$str\n";
	$found = 0;
	$progr = 0;

	if (file_exists("/tmp/daniel_immortal_bicomplex")) system("rm /tmp/daniel_immortal_bicomplex");

	$a = $distance;
	for ($b=0; $b<=$distance; $b++) {
	for ($c=0; $c<=$distance; $c++) {
	for ($d=0; $d<=$distance; $d++) {
	if ($d == 0) echo "\r".round(memory_get_usage()/1024/1024, 2)." MiB ($found found); ($a), $b, $c, $d                                ";
	for ($power=2; $power<=MAX_POWER; $power++) {
	_go($a, $b, $c, $d, $power, $distance, $progr, $found);
	_go($b, $a, $c, $d, $power, $distance, $progr, $found);
	_go($b, $c, $a, $d, $power, $distance, $progr, $found);
	_go($b, $c, $d, $a, $power, $distance, $progr, $found);
	}
	}
	}
	}

#	file_put_contents(__DIR__.'/result.txt', "$str\n", FILE_APPEND);
#	chdir(__DIR__);
#	system("cat /tmp/daniel_immortal_bicomplex | sort | uniq >> result.txt");
#	system("rm /tmp/daniel_immortal_bicomplex");

	system("mv /tmp/daniel_immortal_bicomplex result_distance{$distance}_maxpower".MAX_POWER.".txt");
	system("sort -us -o result_distance{$distance}_maxpower".MAX_POWER."_unique.txt result_distance{$distance}_maxpower".MAX_POWER.".txt");

	file_put_contents(__DIR__.'/savestate.txt', $distance+1);
}