function principal_avalia_hier3_adap(parametros,plota)
%Encontra uma aproximacao de uma funcao dada com base tensorial
%hierarquica adaptativa.

assert(parametros.funcao > 0);
%determina o passo dos niveis
passov = passo_nivel(parametros);

assert(parametros.dim_dominio == 2);
%ajusta o tamanho do dominio conforme o tipo de arvore
area = prod(parametros.tam_dominio);
if parametros.tipoarvore == 1
    LD = [sqrt(sqrt(2)), sqrt(sqrt(0.5))]*sqrt(area);
else
    LD = [1, 1]*sqrt(area);
end
fprintf(2,'tamanho do dominio ajustado para %.3f x %.3f\n',LD(1),LD(2));
parametros.tam_dominio = LD;
CD = 0.5*LD;

%define a malha completa no nivel maximo
fprintf(2,'definindo a malha regular completa do nivel maximo\n');
NC_fin = tamanho_malha_hierarquica(parametros.max_nivel);
nc_fin = prod(NC_fin);
fprintf(2,'numero de celulas no nivel maximo %d x %d = %d\n',...
    NC_fin(1),NC_fin(2),nc_fin);

%calcula o raio de expansao RS do suporte alem da celula central (em celulas) 
rsup_rel = suporte_da_funcao_mae(parametros.tipobase)*parametros.raio;
rs = ceil(rsup_rel - [0.5, 0.5]); 
RS = [rs rs];

%determina a lista p_comp de pontos de amostragem e seus indices PR_comp na
%grade de amostragem completa regular de tamanho NP_comp. No caso de pontos
%regulares, PR_comp eh vazia. 

if parametros.n_pontos == 0
    %determina pontos por celula da malha maxima
    
    KP_fin = calcula_pontos_por_celula(parametros.pts_por_celula,LD,NC_fin);
    fprintf(2,'numero de pontos por celula no nivel maximo %d x %d = %d\n',...
        KP_fin(1),KP_fin(2),KP_fin(1)*KP_fin(2));
    NP_comp = KP_fin .* NC_fin;
    p_comp = constroi_grade_de_amostras(CD,LD,NP_comp);
    np_comp = size(p_comp,1);
    assert(np_comp == prod(NP_comp));
    PR_comp = 1:np_comp;
else
    NP_comp = [NaN NaN];
    p_comp = obtem_pontos_irregulares(parametros.n_pontos,CD,LD);
    np_comp = size(p_comp,1);
    assert(np_comp == parametros.n_pontos)
    PR_comp = zeros(1,0);
end

fprintf(2,'numero de pontos de amostragem  = %d\n',np_comp);

%tamanho da grade de centros varia com o nivel
assert(isnan(parametros.tam_grade_centros(1)))
assert(isnan(parametros.tam_grade_centros(2)))

%tamanho da malha inicial completa
NC_ini = tamanho_malha_hierarquica(parametros.min_nivel);

assert(max(mod(NC_fin,NC_ini)) == 0)

assert(parametros.raios == parametros.raio)
raio_rel = parametros.raio;

%define as pastas para os plotes
sigla_metodo = ['MQ';'SS'; 'SO'; 'WQ'];
pasta_raiz = cria_pasta(sprintf(...
    '~/Desktop/dout2012/Texto/Comum/E4-AH33-F%02d',parametros.funcao));
dir1 = sprintf('arv%d-tb%02d-ten%d-dob%d-gs%d-%s',parametros.tipoarvore,...
    parametros.tipobase,parametros.tensorial,1,...
    0,sigla_metodo(parametros.metodo,:));
if parametros.n_pontos == 0
    dir2 = sprintf('NX%05d-NY%05d-r%08.5f',NP_comp(1),NP_comp(2),raio_rel);
else
    dir2 = sprintf('NP%08d-r%08.5f',parametros.n_pontos,raio_rel);
end
pasta_global = cria_pasta(sprintf('%s/%s/%s',...
    pasta_raiz,dir1,dir2));

%constroi a matriz f_comp da funcao a aproximar
%f_comp(i,j) valor da funcao numero j ( = 1) no ponto de amostragem i
f_comp = discretiza_funcao(parametros,p_comp,CD,LD);
assert(size(f_comp,2) == 1);
nf = 1;
if plota.funcao
    figure(9);
    for j = 1:nf
        plota_funcao(p_comp,f_comp(:,j),true,CD,LD,NP_comp,NaN,NaN);
        %define_zlim2(alpha,f,s);
        title(sprintf('funcao'),'FontSize',20);
        hold off
        if plota.salva_eps
            saveas(gcf,sprintf('%s/f%04d.eps',pasta_raiz,j));
        end
        %         pause(1)
    end
end

arq_global = fopen(sprintf('%s/parametros.tex',pasta_global),'wt');
fprintf(arq_global,'\\newcommand{\\PARMtamdomX}{%d}\n',parametros.tam_dominio(1));
fprintf(arq_global,'\\newcommand{\\PARMtamdomY}{%d}\n',parametros.tam_dominio(2));
fprintf(arq_global,'\\newcommand{\\PARMtensorial}{%d}\n',parametros.tensorial);
%fprintf(arq_global,'\\newcommand{\\PARMdobrado}{%d}\n',parametros.dobrado);
if parametros.n_pontos == 0
    fprintf(arq_global,'\\newcommand{\\PARMamostrX}{%d}\n',NP_comp(1));
    fprintf(arq_global,'\\newcommand{\\PARMamostrY}{%d}\n',NP_comp(2));
end
fprintf(arq_global,'\\newcommand{\\PARMamostras}{%d}\n',np_comp);
fprintf(arq_global,'\\newcommand{\\PARMraio}{%s}\n',...
    converte_raio_tex(raio_rel));

NC_ant = [NaN NaN];

%base_h eh a base hierarquica (multinivel) atual
%nb_h eh numero total de elementos na base hierarquica
%alpha_h eh o vetor de coeficientes para a base base_h
%A_h eh a matriz desta base calculada nos pontos p_comp
%inib_h(k) eh o comecco do nivel k em alpha_h, base_h
%tamb_h(k) eh o numero de elementos da base no nivel k
%malha_h(k) eh a malha com as celulas ativas do nivel k
profile on
nb_h = 0;
alpha_h = zeros(nb_h,nf);
A_h = sparse(np_comp,0);
beta = 1; %reducao de tolerancia por nivel 
theta = 1; %tolerancia relativa para a escolha das celulas criticas
for iv = parametros.min_nivel:passov:parametros.max_nivel
    %disp('estamos no pause')
    %pause
    fprintf(2,'=== nivel %d ======================================\n',iv);
    
    %calcula a tolerancia a usar neste nivel
    tol = parametros.tol*beta^(parametros.max_nivel - iv);
    fprintf(2,'tolerancia = %8.6f \n',tol);
    
    plota.janela_base = 20;
    pasta_nivel = cria_pasta(sprintf('%s/iv%02d',pasta_global,iv));
    
    %constroi a malha do nivel iv (malha_v)
    if iv == parametros.min_nivel
        %constroi a malha inicial
        NC_v = NC_ini;
        malha_v = constroi_malha_para_pontos(CD,LD,NC_v,p_comp,RS);
        e_comp = f_comp; 
        s_comp = zeros(np_comp,nf);
    else
        %refina a malha anterior, que deve ter todas as celulas ativas (que
        %podem ainda ser modificadas).
        KF = numero_de_filhos(parametros.tipoarvore,parametros.dim_dominio,...
            iv-passov);
        fprintf(2,'numero de filhos por celula da malha anterior = [%d %d] \n',...
            KF(1),KF(2));
        NC_v = NC_ant .* KF;
        malha_v = refina_malha(malha_ant,KF,p_comp);
    end
    fprintf(2,'numero de celulas na malha completa %d x %d = %d\n',...
        NC_v(1),NC_v(2),prod(NC_v));
    
    estagio = sprintf('nivel %02d (preb)',iv);
    
    %determina celulas onde hah erro excessivo  
    K_v = determina_celulas_criticas(malha_v,p_comp,e_comp,theta*tol);
    nk_v = size(K_v,1);
    fprintf(2,'numero de celulas criticas = %d\n',nk_v);
    
    %calcula o conjunto de celulas dos suportes da base
    S_v = expande_conjunto_de_celulas(K_v,NC_v,RS);
    ns_v = size(S_v,1);
    fprintf('numero de celulas no suporte = %d\n',ns_v);
    
    %calcula o conjunto de celulas ativas que podem se afetadas daqui para
    %frente
    U_v = expande_conjunto_de_celulas(S_v,NC_v,RS);
    nu_v = size(U_v,1);
    fprintf('numero de celulas ativas = %d\n',nu_v);
    
    figure(plota.janela_base + 15)
    plota_celulas_KSU(K_v,S_v,U_v,'pre',CD,LD,NC_v,estagio,plota,pasta_nivel)
    
    %elimina celulas que ficaram inativas
    malha_v = encolhe_malha(malha_v,U_v);
    
    %nc_v = size(malha_v.celulas,2);
    nc_v = nnz(malha_v.ixc);
     fprintf(2,'numero de celulas na malha corrente = %d\n',nc_v);
    assert(nc_v == nu_v)
    
    if nc_v == 0
        %fprintf(arq_global,'\\newcommand{\\PARMnbh}{%d}\n',nb_h);
        break
    end
   
    raio_abs_v = (LD./NC_v)*raio_rel;
    
    fprintf(2,'raio absoluto do nivel %d = [%8.5f %8.5f]\n',...
        iv,raio_abs_v(1),raio_abs_v(2));
    
    %constroi a prebase do nivel iv
    base_v = constroi_base(malha_v,K_v,raio_abs_v,parametros.tipobase);
    nb_v = size(base_v.elementos,2);
    fprintf(2,'tamanho da pre-base deste nivel %d\n',nb_v);
    assert(nb_v == nk_v)
    
   
    
    %obtem indices PR_sup dos pontos de amostragem das celulas do suporte,
    %no vetor p_comp de todos os pontos de amostragem
    PR_sup = obtem_indices_dos_pontos_das_celulas(malha_v,S_v);
    p_sup = p_comp(PR_sup,:);
    np_sup = size(p_sup,1);
    e_sup = e_comp(PR_sup,:);
    fprintf(2,'numero de pontos no suporte da prebase = %d\n', np_sup);
    
    
    %obtem a matriz A_sup da prebase base_v calculada nos pontos p_sup
    
    [A0_comp,A_sup] = obtem_matriz_A(parametros,p_sup,malha_v,base_v,...
        parametros.tensorial,NP_comp,np_comp,PR_sup);
    assert(size(A_sup,1) == np_sup)
    assert(size(A_sup,2) == nb_v)
    
    if parametros.metodo == 2
        alpha_v = calcula_aprox_inicial_IC_hier3(parametros,...
            malha_v,base_v,p_comp,e_comp,plota,estagio,pasta_nivel);
    else
        assert((parametros.metodo == 1)|| (parametros.metodo == 3))
        alpha_v = calcula_aprox_inicial_MQ_hier3(parametros,...
            malha_v,base_v,A0_comp,A_sup,NP_comp,p_comp,e_comp,...
            PR_sup,plota,estagio,pasta_nivel);
    end
    
    arq = fopen(sprintf('%s/parametros.tex',pasta_nivel),'wt');
    fprintf(arq,'\\newcommand{\\PARMraioX}{%s}\n',...
        converte_raio_tex(raio_abs_v(1)));
    fprintf(arq,'\\newcommand{\\PARMraioY}{%s}\n',...
        converte_raio_tex(raio_abs_v(2)));
    fprintf(arq,'\\newcommand{\\PARMcentrosX}{%d}\n',NC_v(1));
    fprintf(arq,'\\newcommand{\\PARMcentrosY}{%d}\n',NC_v(2));
    fprintf(arq,'\\newcommand{\\PARMdimAnivel}{%d}\n',nb_v);
    
    es_sup = A_sup*alpha_v;
    ee_sup = e_sup - es_sup;
    erro_sup_L2 = norma_erro_L2(ee_sup,np_sup);
    erro_sup_Linf = norma_erro_Linf(ee_sup);
    fprintf(2,'erro com prebase %9.6f (L2) %9.6f (Linf)\n',...
        erro_sup_L2, erro_sup_Linf);
    fprintf(arq,'\\newcommand{\\PARMerroHilb}{%s}\n',...
        converte_erro_tex(erro_sup_L2));
    fprintf(arq,'\\newcommand{\\PARMerroLinf}{%s}\n',...
        converte_erro_tex(erro_sup_Linf));
    
    if parametros.adap_por_nivel 
        fprintf(2,'--- aproximacao adaptativa no nivel %d -------------\n',iv);
        estagio = sprintf('nivel %02d (adap)',iv);
        plota.janela_base = 20;  %reusa as janelas
        pasta_nivel = cria_pasta(sprintf('%s-a',...
            pasta_nivel));
        plota.base_dual = false;
        
        %calcula a base reduzida base_a, A_a,alpha_a, etc a partir da base 
        %base_v, tudo  calculado nos pontos p_sup.
        [ix_a,base_a,A_a,alpha_a] = ...
            escolhe_sub_base(parametros,base_v,A_sup,alpha_v,ee_sup);
        assert(size(A_a,1) == np_sup)
        nb_a = size(A_a,2);
        fprintf(2,'tam base adaptativa do nivel = %d\n',nb_a);
       
        
                
        
        %calcula a aproximacao es_a e o residuo ee_a com coeficientes 
        %originais e base reduzida, nos pontos p_sup
        es_a = A_a*alpha_a;
        if iv == 17
            size(es_a)
            size(e_sup)
        end
        ee_a = e_sup - es_a;
        erro_a_L2 = norma_erro_L2(ee_a,np_sup);
        erro_a_Linf = norma_erro_Linf(ee_a);
        fprintf(2,'erro com base reduzida e alphas originais %9.6f (L2) %9.6f (Linf)\n',...
            erro_a_L2, erro_a_Linf);
        arq_a = fopen(sprintf('%s/paramadapt.tex',pasta_nivel),'wt');
        fprintf(arq_a,'\\newcommand{\\PARMerroaHilb}{%s}\n',...
            converte_erro_tex(erro_a_L2));
        fprintf(arq_a,'\\newcommand{\\PARMerroaLinf}{%s}\n',...
            converte_erro_tex(erro_a_Linf));
        fprintf(arq_a,'\\newcommand{\\PARMtol}{%s}\n',...
            converte_erro_tex(tol));
        fprintf(arq_a,'\\newcommand{\\PARMnba}{%d}\n',nb_a);
        fclose(arq_a);
        
        %substitui a prebase pela base reduzida
        base_v = base_a;
        nb_v = nb_a;
        A_sup = A_a;
        alpha_v = alpha_a;
        es_sup = es_a;
        ee_sup = ee_a;
        erro_sup_L2 = erro_a_L2;
        erro_sup_Linf = erro_a_Linf;
        
    end
    fclose(arq);
    assert(size(base_v.elementos,2) == nb_v)
    
    %insere a aproximacao e o residuo nos vetores globais s_comp, e_comp
    s_comp(PR_sup,:) = s_comp(PR_sup,:) + es_sup;
    e_comp(PR_sup,:) = ee_sup;
    
    %plota a aproximacao acumulada incluindo este nivel
    if plota.aproximacao 
        assert(nf == 1)
        estagio = sprintf('nivel %02d (acum)',iv);
        nomef = 'funcao';
        %por enquanto imprime na mesma janela e na pasta global
        janela_base = 40;
        pasta = pasta_global;
        figure(5+janela_base);
        plota_funcao(p_comp,s_comp,true,CD,LD,NP_comp,NaN,NaN);
        title(sprintf('aproximacao da %s %s',...
            nomef,estagio),'FontSize',20);
        hold off
        if plota.salva_eps
            saveas(gcf,sprintf('%s/s_v.eps',pasta));
        end
        
        figure(6+janela_base);
        if plota.erro_curva_nivel
            plot_contour_funcao(parametros,p_comp,e_comp,CD,LD,NP_comp);
        else
            emax = max([plota.emax, max(abs(e_comp))]);
            plota_funcao(p_comp,e_comp,false,CD,LD,NP_comp,-emax,emax);
        end
        title(sprintf('erro da %s %s',...
            nomef,estagio),'FontSize',20);
        hold off
        if plota.salva_eps
            saveas(gcf,sprintf('%s/e_comp.eps',pasta));
        end
    end
    
    estagio = sprintf('nivel %02d (final)',iv);
    
    %recalcula o conjunto de celulas da base do nivel
    K_b = obtem_posicoes_dos_elementos(base_v);
    nk_b = size(K_b,1);
    fprintf(2,'numero final de celulas da base do nivel = %d\n',nk_b);
    assert(nk_b == nb_v)
    
    %calcula um conjunto de celulas que contem todas as celulas criticas 
    %do nivel seguinte
    S_b = expande_conjunto_de_celulas(K_b,NC_v,RS);
    ns_b = size(S_b,1);
    %acrescenta as celulas criticas iniciais
    S_b(ns_b+1:ns_b+nk_v,:) = K_v;
    S_b = posicoes_unicas(S_b);
    ns_b = size(S_b,1);
    fprintf('numero final de celulas no "suporte do nivel = %d\n',ns_b); 
    
    %calcula o conjunto das celulas onde a aproximacao ainda pode ser mudada
    %nos niveis seguintes
    U_b = expande_conjunto_de_celulas(S_b,NC_v,RS);
    nu_b = size(U_b,1);
    fprintf('numero final de celulas ativas do nivel = %d\n',nu_b);
    
    figure(plota.janela_base + 16)
    plota_celulas_KSU(K_b,S_b,U_b,'fim',CD,LD,NC_v,estagio,plota,pasta_nivel)
    
    malha_s = encolhe_malha(malha_v,U_b);
    
    %acrescenta elementos do nivel iv na base multinivel
    malha_h(iv) = malha_s;
     ini_v = nb_h + 1;
     fim_v = ini_v + nb_v - 1;
     A_h(PR_sup, ini_v:fim_v) = A_sup;
    if nb_h == 0
        base_h = base_v;
        alpha_h = alpha_v;
        inib_h(iv) = 1;
    else 
        base_h.elementos(1,ini_v:fim_v) = base_v.elementos(1,:);
        %alpha_h(1,ini_v:fim_v) = alpha_v(1,:);
         alpha_h(ini_v:fim_v,1) = alpha_v(:,1);
        base_h.uniforme = false;
        inib_h(iv) = ini_v;
    end
    tamb_h(iv) = nb_v;
    nb_h = nb_h + nb_v; 
    fprintf(2,'tamanho da base hierarquica acumulada = %d\n', nb_h);
    
    figure(100); clf;
                plota_centros_hier2(base_h,p_sup,A_h,malha_h,NP_comp,...
                    parametros,parametros.min_nivel,passov,iv,inib_h,tamb_h,...
                    pasta_nivel,plota,estagio)
    
    %calcula o conjunto de celulas ativas (que podem se afetadas daqui para
    %frente)
    U_b = expande_conjunto_de_celulas(S_b,NC_v,RS);
    nu_b = size(U_b,1);
    fprintf('numero final de celulas ativas = %d\n',nu_b);
    
    %encolhe a malha_v para ter apenas as celulas ainda ativas.
    malha_ant = encolhe_malha(malha_v,U_b);
    NC_ant = NC_v;
end
%recalcula erro da aproximacao acumulada
erro_comp_L2 = norma_erro_L2(e_comp,np_comp);
erro_comp_Linf = norma_erro_Linf(e_comp);
normae = norm(e_comp,inf);
fprintf(2,'erro com base acumulada %9.6f (L2) %9.6f (Linf)\n',...
    erro_comp_L2, erro_comp_Linf);

fprintf(arq_global,'\\newcommand{\\PARMerrocompHilb}{%s}\n',...
    converte_erro_tex(erro_comp_L2));
fprintf(arq_global,'\\newcommand{\\PARMerrocompLinf}{%s}\n',...
    converte_erro_tex(erro_comp_Linf));

fprintf(arq_global,'\\newcommand{\\PARMnormaeLinf}{%s}\n',...
    converte_erro_tex(normae));
fprintf(arq_global,'\\newcommand{\\PARMnbh}{%d}\n',nb_h);
fclose(arq_global);
profile viewer
%deveria tentar achar a aproximacao e/ou plotar com a base inteira