%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% optimisation topologique
% en contraintes planes
% version + elaboree
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% initialisation :
% c = tenseur des contraintes
% hm <= h = densite variable <= hp
% reduc = proportion retiree a chaque
% iteration
% a prendre .85 <= reduc <= .98
% propfin = proportion finale matiere
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear all; global hm hp b p e t a f
load pont; reduc = .9; propfin = .4;
[aux,nbt] =  size(t);
cplein = fairec(c,nbt); c = cplein;
h = ones(1,nbt); hm = 1e-3; hp = 1;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% suppression des elements
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

plein = 1; figure(1); k = 1;

while plein >= propfin
 plein = plein*reduc; vide = 1-plein;
 disp('enleve, plein='); disp(plein);
 [c,h] = retire(c,cplein,h,vide);

 figure(2);
 pdeplot(p,e,t,'xydata',h,...
 'xystyle','flat','mesh','off');
 pause(1);
 %impr=['print -deps pont2' int2str(k)];
 %eval(impr); k = k+1;
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%iteration a volume constant
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

for k = 1:0
 disp(k);
 [c,h] = retire(c,cplein,h,vide);
 figure(2);
 pdeplot(p,e,t,'xydata',h,...
 'xystyle','flat','mesh','off');
 pause(.1);
end

figure(1);
voirmis;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
