function	my_moo

%Calculates multiplexor approximation 
%proposed in quant-ph/0412072

NB=4;   %NB=positive integer, number of bits
example =1;

if (NB<1)
	error("NB is less than 1");
end
len_phi = 2^(NB-1);
phi=zeros(len_phi, 1);
switch example
	case (1)  %phi in increasing order
		%rand("seed", 1.27);
		phi = sort(rand(len_phi,1));
	case (2)  %phi in decreasing order
		%rand("seed", 1.27);
		phi = flipud(sort(rand(len_phi,1)));
	case (3)
		if(len_phi!=8)
			error("this example requires NB=4");
		end
		phi=[.31; .31;.3100001; .31005; .5;.5;.5; .5];
	case (4)
		if(len_phi!=16)
			error("this example requires NB=5");
		end
		phi=[0;0;0;0;0;0;0;1;1;1;1;1;1;1;1;1];
	otherwise
		error("example number is out of range");
end
max_phi =norm(phi, Inf); % all phi components are non-negative

err_file = fopen ("out_error.txt", "w", "native");
phi_file = fopen ("out_phis.txt", "w", "native");

for i=1:len_phi
	fprintf(phi_file, "phi(%d)=\t%11.9f\n",i, phi(i));
end

fprintf(err_file, "error as function of (permutation\\delta_B)\n");
for del_B=0:(NB-1 )
	fprintf(err_file, "\t%9d", del_B);
end
fprintf(err_file, "\n");

more = false;
pi_B = (1: NB-1);
perm_num=0;
while (1)
	 [ pi_B_new, more_new ] = perm_lex_next ( NB-1, pi_B, more );
	if (more_new)
		pi_B =  pi_B_new;
		more = true;
		perm_num++;

		fprintf(phi_file, "-----------------------\n");
		fprintf(phi_file, "permutation %d = (", perm_num);
		for   i = 1 : NB-2
			fprintf (phi_file, "%d,", pi_B(i) );
		end
		fprintf (phi_file, "%d)\n", pi_B(NB-1));
		fprintf(phi_file, "delta_B, phi_prime=\n");

		fprintf(err_file, "%4d", perm_num);

		for del_B=0:(NB-1 )
			phi_pr = approx_phi(phi, pi_B, del_B, NB);

			fprintf(phi_file, "%d", del_B);
			for i=1:len_phi
				comp= phi_pr(i);
				fprintf(phi_file,"\t%5.3f", comp);
			end
			fprintf(phi_file, "\n");

			err= norm(phi - phi_pr, Inf);
			fprintf(err_file, "\t%.3e", err);

		end
		fprintf(err_file, "\n");
	else
		break;
	end
end
fclose(err_file);
fclose(phi_file);