read "itall.txt":

MAKEBASEEIGENS:=0:
MAKEALLEIGENS:=0:
SAVENORMS:=0:
NORMALIZEEIGENS:=0:
TESTSTUFF:=0:
	CHECKORTHONORMAL:=1:
	TOLNORMALITY:=.0000000001:
	TOLORTHOGONALITY:=.0000000001:
	CHECKEIGENVALS:=1:
	TOLEIGENVALS:=.0000000001:
FATTENEIGENS:=0:
MAKEPSIS:=0:
	PSIMAKINGTOLERANCE:=0;
MAKEDIRCHLET:=1:
	BMTOL:=.000000000000001:
CHECKPHATEIGENVALS:=1:
	

lvl:=3;
phatlvl:=8;
path:=cat("../data/level",lvl,"/");
inpath:=cat("../data/level",lvl,"/");
outpath:=cat("../data/level",lvl,"/");
numEigens:=(3^(lvl+1)-3)/2;

T:=gasket(phatlvl):


if(MAKEBASEEIGENS=1) then
	makeOrthoEigen5(T,lvl,outpath);
	makeOrthoEigen6(T,lvl,outpath);
	print("All done making basic orthonormal Eigen Set");
end if:



if(MAKEALLEIGENS=1) then
	genalllevel(T,lvl,path);
	print("All done making full orthogonal Eigen Set");
end if:



if (SAVENORMS=1) then
	fd:=fopen(cat(path,lvl,"squarednorms.txt"),WRITE);

	for i from 1 to numEigens do
        	readFunc(T,12,cat(path,lvl,"eigen.",i));
        	fprintf(fd,"%20.20g\n",discInnerProduct(T,lvl,12,12));
	end do;

	fclose(fd);
	print("All done Saving Norms");
end if:



if (NORMALIZEEIGENS=1) then
	for i from 1 to numEigens do
       		readFunc(T,12,cat(path,lvl,"eigen.",i));
	        normalizeFunc(T,lvl,12);
        	saveFunc(T,12,lvl,cat(path,lvl,"eigen.",i));
	end do;
	print("All done Normalizing");
end if:



if (TESTSTUFF=1) then
	startEigen:=12;
	tol2:=TOLEIGENVALS;
	tol1:=TOLORTHOGONALITY;
	tol0:=TOLNORMALITY;

	for i from 1 to numEigens do
	        readFunc(T,i+startEigen-1,cat(path,lvl,"eigen.",i));
	end do;

	if (CHECKORTHONORMAL=1) then
		for i from 1 to numEigens do
		        for j from i to numEigens do
        		        m:=discInnerProduct(T,lvl,i+startEigen-1,j+startEigen-1);
                        
		                if (i=j) then
		                        if(abs(m-1) > tol0) then
		                        printf("(%d,%d)=%20.20g\n",i,j,m);
        	        	        end if;
        		        else    if(abs(m) > tol1) then                
                	        	printf("(%d,%d)=%20.20g\n",i,j,m);
                		        end if;
        		        end if;
	
       			 end do;
		end do;
		print("All done Checking Orthonormality");	
	end if;

	if(CHECKEIGENVALS=1) then

		fd:=fopen(cat(path,lvl,"lambdalist.txt"), READ);
		for i from 1 to numEigens do
        		discLaplac(T,lvl,i+startEigen-1,startEigen+numEigens);
       			funcQuotient(T,lvl,startEigen+numEigens,i+startEigen-1,startEigen+numEigens);
        		lambda:=fscanf(fd,"%g\n")[1];
        		lambda:=-1*lambda;
       			m:=compareFuncList(T,lvl,startEigen+numEigens,[0,lambda]);
        		if (abs(m) > tol2) then 
                		printf("max variance of (laplacian(%d)/%d) =%20.20g\n",i,i,m);
        		end if;
		end do;

		fclose(fd);

		print("All done Checking Eigens");

	end if;

end if:



if (FATTENEIGENS=1) then

	fd:=fopen(cat(inpath,lvl,"lambdalist.txt"), READ);
	df:=fopen(cat(outpath,phatlvl,"phat",lvl,"lambdalist.txt"), WRITE);
	
	for i from 1 to numEigens do
	        lambda:=fscanf(fd,"%g\n")[1];
	        readFunc(T,12,cat(inpath,lvl,"eigen.",i));
	        if (lambda=6) then 
        	        if (phatlvl > lvl) then
        	                fixedeps:=[1,op(makeEpsilons(phatlvl - lvl -1))];
        	                newAlgoWrapper(T,lvl,12,lambda,fixedeps);
        	        	newlambda:=makeLambda(lambda,fixedeps);
			end if;
       		else    newAlgoWrapper(T,lvl,12,lambda,makeEpsilons(phatlvl - lvl));
			newlambda:=makeLambda(lambda,makeEpsilons(phatlvl - lvl));
        	end if;
	        fprintf(df,"%20.20g\n",newlambda);
		saveFunc(T,12,phatlvl,cat(outpath,phatlvl,"phat",lvl,"eigen.",i));
	end do;
	
	fclose(df);
	fclose(fd);
	print("All done Fattening Eigens");	
	
end if:



if (MAKEPSIS=1) then
	tol:=PSIMAKINGTOLERANCE;

	vl:=ourVertices(lvl):
	vl:=subsop(1=NULL,2=NULL,3=NULL,vl):
	vp:=ourVertices(phatlvl):

	for i from 1 to nops(vl) do
	        clearFunc(T,phatlvl,13);
	        for j from 1 to nops(vl) do

			readFunc(T,12,cat(inpath,phatlvl,"phat",lvl,"eigen.",j));
	                co:=T[vl[i],12];
	                if (abs(co) > tol) then
        	                funcSumWithConstant(T,phatlvl,co,13,12,13);
        	        end if;
       		end do;
        
		saveFunc(T,13,phatlvl,cat(outpath,phatlvl,"phat",lvl,"psi.",i));
	end do;

	print("All done Making Psis");
end if:


if (MAKEDIRCHLET=1) then

	vl:=ourVertices(lvl):
	vl:=subsop(1=NULL,2=NULL,3=NULL,vl):
	vp:=ourVertices(phatlvl):

	for i from 1 to nops(vl) do
	        clearFunc(T,phatlvl,13);
		fd:=fopen(cat(path,lvl,"lambdalist.txt"), READ);
	        for j from 1 to nops(vl) do

			ulambda:=fscanf(fd,"%g\n")[1];

			readFunc(T,12,cat(inpath,phatlvl,"phat",lvl,"eigen.",j));
	                co:=T[vl[i],12]*(3^(lvl+1))/(2*calcbm(ulambda,[],BMTOL));
	                if (abs(co) > 0) then
        	                funcSumWithConstant(T,phatlvl,co,13,12,13);
        	        end if;
       		end do;
		fclose(fd);
        
		saveFunc(T,13,phatlvl,cat(outpath,phatlvl,"phat",lvl,"dirichlet.",i));
	end do;

	print("All done Making Psis");
end if:


if(CHECKPHATEIGENVALS=1) then

	df:=fopen(cat(path,phatlvl,"phat",lvl,"lambdaerrs.txt"),WRITE);
	fd:=fopen(cat(path,phatlvl,"phat",lvl,"lambdalist.txt"), READ);
	
	for i from 1 to numEigens do
	   
		readFunc(T,startEigen,cat(path,phatlvl,"phat",lvl,"eigen.",i));
	
	       	discLaplac(T,phatlvl,startEigen,startEigen+1);
	        funcQuotient(T,phatlvl,startEigen+1,startEigen,startEigen+2);
	        lambda:=fscanf(fd,"%g\n")[1];
	        lambda:=-1*lambda;
	        m:=compareFuncList(T,phatlvl,startEigen+2,[0,lambda]);
	        fprintf(df,"Laplacian error for function %d =\t%20.20g\n",i,m);
	                
	end do;

	fclose(fd);
	fclose(df);

	print("All done Checking Phat Eigens");

end if;


print("Have a nice Day");