// Copyright (C) 2004-2011 - J.-Ph. Chancelier. // // This library is free software; you can redistribute it and/or // modify it under the terms of the GNU General Public // License as published by the Free Software Foundation; either // version 2 of the License, or (at your option) any later version. // // This library is distributed in the hope that it will be useful, // but WITHOUT ANY WARRANTY; without even the implied warranty of // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU // General Public License for more details. // // You should have received a copy of the GNU General Public // License along with this library; if not, write to the // Free Software Foundation, Inc., 59 Temple Place - Suite 330, // Boston, MA 02111-1307, USA. // // Source: "Introduction à Scilab", Chancelier, 2004 // Chapitre 2: Programmer : Arbres et matrices de ramification // // With adaptations for Scilab v5, 2011, Michael Baudin stacksize("max"); // 2.0.5 Matrices de ramifications function [N,mat]=matrice_ramification(N,val) select val case 1 then mat=eye(N,N); case 2 then mat= (1/2).^[1:N-1]; mat=ones(N-1,1)*mat ; mat=tril(mat); d=diag(diag(mat),1);d($,:)=[]; mat = [mat,zeros(N-1,1)]+d; mat = [eye(1,N);mat]; end endfunction [N,mat]=matrice_ramification(3,2); mat // 2.0.6 Tirage aleatoire des ordres des aretes filles function [ordreg,ordred]=calcul_ordre(k) o1=k; o2=grand(1,'markov',mat(k,1:k),1) if o2==k then o1=k-1; o2= k-1; end ale= grand(1,'uin',1,2) ; if ale==1, ordreg=o1; ordred=o2; else ordreg=o2; ordred=o1; end endfunction function [L]=construire(k) if k==1 then L=list(1); else [ordreg,ordred]=calcul_ordre(k); L=list(k,construire(ordreg),construire(ordred)); end endfunction N=2; [N,mat]=matrice_ramification(N,1); L=construire(N) // 2.0.8 Utilisation des graphes Scilab pour visualiser l'arbre function [L,nmax]=construire(k,numero,prof) if k==1 then L=list([1,numero,prof]); nmax=numero; return; end [ordreg,ordred]=calcul_ordre(k); [Lg,nmg]=construire(ordreg,numero,prof+1) [Ld,nmd]=construire(ordred,nmg+2,prof+1); L=list([k,nmg+1,prof],Lg,Ld); nmax=nmd; endfunction function [M,pos]=construire_matrice(L,M,pos) if size(L)== 3 then M=[M,[L(1)(2);L(2)(1)(2)],[L(1)(2);L(3)(1)(2)]]; [M,pos]=construire_matrice(L(2),M,pos); [M,pos]=construire_matrice(L(3),M,pos); end pos=[pos,L(1)'] endfunction N=4;val=1; [N,mat]=matrice_ramification(N,val); [L,nmax]=construire(N,1,1); if ( %f ) then // TODO : fixme [g,rect]= construire_graphe(L,nmax); rep=[2 1 1 2 1 2 2 2 2 2 2 2 2]; plot_graph(g,rep,rect); end // 2.1 Representation graphique directe function []=empiler(x) // Rajoute un vecteur ligne au sommet de la pile, la taille de x doit ^etre fixe global %pile; n=length(x);buf_size=1000; // buf_size nous donne la taille des blcos if typeof(%pile)=='constant' then // Creation de la pile %pile=tlist(['pile','val','pos','libres','inc'],0*ones(buf_size,n),1,buf_size,buf_size); end if %pile.libres > 1 then %pile.val(%pile.pos ,:)=x; // On stocke x %pile.pos = %pile.pos + 1; // On incremente la position libre %pile.libres = %pile.libres - 1; // On decremente le nombre de positions libre else %pile.val=[%pile.val;0*ones(%pile.inc,n)]; // On rajoute un block dans la pile %pile.libres = %pile.libres + %pile.inc; // On mets a jours le nombre de positions libres empiler(x); // On stocke x end endfunction function m=valeurs() // On extrait de la pile la matrice des valeurs stockees global %pile; m=%pile.val(1:%pile.pos-1,:); // La derniere valeur stockee est sur la ligne pos-1 endfunction function []=calculer(cp,cdd,k) // calcul les ar^etes fille du segment cdd 7! cp d'ordre k empiler([cdd(1),cp(1),cdd(2),cp(2),k]); // On empile les informations sur le segment if k==1,return;end // Arr^et sur une feuille [ordreg,ordred]=calcul_ordre(k); // Calcul des ordres des deux branches filles l1=l(ordreg),l2=l(ordred); // Les longueurs des branches filles delta1=del(ordreg)*%pi/180,delta2=del(ordred)*%pi/180; // les angles de rotation if ordreg==k then delta1=delta1/5; // la branche principale tourne moins end if ordred==k then delta2=delta2/5; // la branche principale tourne moins end [cp1,cp2]=newp(cp,cdd) // calcul des deux nouveaux noeuds calculer(cp1,cp,ordreg); // calcul du nouveau segment qui suit cp 7! cd1 calculer(cp2,cp,ordred); // calcul du nouveau segment qui suit cp 7! cd2 endfunction function [cp1,cp2]=newp(cp,cdd) // Rotations pour trouver les coordonnees des deux noeud fils ddir=cp-cdd; ddir = ddir/norm(ddir,2) diro=[ddir(2),-ddir(1)] cp1= cp + l1*cos(delta1)*ddir + l1*sin(delta1)*diro cp2= cp + l2*cos(delta2)*ddir - l2*sin(delta2)*diro endfunction clearglobal %pile; N=10; val=1; [n,mat]=matrice_ramification(N,val); l=[1:N]/N; // les longueurs a utiliser pour les aretes de nombre i del=20*ones(1,N); // l'angle a utiliser pour les aretes de nombre i scf(); nc=N; cmap=[((1:nc).^2)/nc;1:nc;((1:nc).^2)/nc]'/nc; xset("colormap",cmap) // l'angle a utiliser pour les aretes de nombre i calculer([0,0],[0,-1],N); m=valeurs(); x=m(:,1:2); y=m(:,3:4); plot2d(0,0,1,"010"," ",[min(x),min(y),max(x),max(y)]) for i=1:N I=find(m(:,5)==i); xset("thickness",i) // l'epaisseur du trait depend de l'ordre de l'arete. xsegs(x(I,:)',y(I,:)',N-i); end