#looks for smallest determinants of the Psi matrix for irregular sampling
#vertices are shifted to the various neighbors (of a given level) of the vertex they should b

read "itall.txt":

with(LinearAlgebra);

mindet:=proc(T,lvl,funclist,listverts,tol)
	local l, m, n, d, oldlow, cycles, lowdet,i,j,k, newlist, tmpverts;
	
	make_neighbors(T,lvl,T[lvl,-1]);
	l:=nops(listverts);
	
	newverts:=listverts;
	m:=Matrix([seq([seq(T[newverts[i],funclist[j]],j=1..l)],i=1..l)]);
	lowdet:=abs(Determinant(m));
	
	cycles:=0;
	oldlow:=lowdet+10*tol;
	
	while (lowdet < oldlow - tol) do	
		oldlow:=lowdet;
		for i from 1 to l do
			for j from 0 to 4 do
				tmpverts:=subsop(i=T[listverts[i],j],newverts);
				m:=Matrix([seq([seq(T[tmpverts[n],funclist[k]],k=1..l)],n=1..l)]);	
				d:=abs(Determinant(m));
				if (d<lowdet) then
					newverts:=tmpverts;
					lowdet:=d;
				end if;
			end do;
		end do;
		cycles:=cycles+1;
	end do;

	return([lowdet,newverts]);

end;

PRINTVERTS:=true:
PRINTFIRSTLEVEL:=true:

lvl:=4;
phatlvl:=8;
T:=verGasket(phatlvl);
tol:=.000001;
startfunc:=12;

path:=cat("../data/level",lvl,"/");
lverts:=subsop(1=NULL,2=NULL,3=NULL,T[lvl,-1]);

for i from 1 to nops(lverts)  do
	readFunc(T,startfunc+i-1,cat(path,phatlvl,"phat",lvl,"psi.",i));
end do;

fd:=fopen(cat(path,lvl,"mindets.txt"),WRITE);
fprintf(fd,"level\t(base level %d)\tdeterminant\n",lvl);

if (PRINTFIRSTLEVEL) then
	st:=lvl+1:
else	st:=lvl+2:
end if;

for level from st to phatlvl do
	d:=mindet(T,level,[seq(12+i-1,i=1..nops(lverts))],lverts,tol);
	print(cat("mindet at level ",level," is ",d[1]));
	if (PRINTVERTS) then print(cat("at vertices ",d[2])); end if;
	fprintf(fd,"%d\t\t%20.20g\n",level,d[1]);
	if (PRINTVERTS) then fprintf(fd,"%a\n",d[2]); end if;
end do;
fclose(fd);