function [un,p,pmax,untot]=troisnplusun(Ntot,graphe,nmax)

% Valeurs et temps de vol de la suite de Syracuse. 
%
% *********************************************************
% [un,p,pmax,untot]=troisnplusun(N,graphe,nmax) : étude de la suite
% de Syracuse, définie par u_0=N entier non nul et 
% u_{n+1}=u_n/2 si u_n est paire
% u_{n+1}=3u_n+1 sinon
% voir Conjecture de Collatz ou de Syracyuse ou
% https://fr.wikipedia.org/wiki/Conjecture_de_Syracuse
%
%  * Variable d'entrées
%     * Ntot : la valeur maximale de N (qui varie entre 1 et Ntot)
%  * Variables d'entrées optionnelles
%     * graphe : entier dans {0,1}, valant 1 par défaut. 
%       si égal à 1, tracé de graphe temps de vol en fonction de N dans
%       {1,...,Ntot}.
%     * nmax : nombre maximal de calcul, égal à 1e7 par défaut.
%  * Variable de sortie
%     * un : tableau d'entier de taille N, où un(i) est la valeur finale
%     de la suite u_n correspondant à u_0=un(i), en principe égale à 1.
%     * p : tableau d'entier de taille N, où p(i) est le temps de vol de la
%     suite de Syracuse N
%    * pmax : max de p.
%  * Variable de sortie optionnelle
%    * ntot : cellule de taille Ntot donnant toutes les valeurs des suites
%    de Syracuse.
%
% Exemples :
% [un,p,pmax]=troisnplusun(1000);pmax
% [un,p,pmax]=troisnplusun(100000);pmax
%
%
% ************ Fonctions auxiliaires utilisées ************
%
%    aucune
%
% *********************************************************
%
% 2018 by          Jérôme BASTIEN
%                  Université Claude Bernard Lyon I, Polytech, Laboratoire LIBM, Villeurbanne
%                  E-Mail : jerome.bastien@univ-lyon1.fr

if nargin<=1||isempty(graphe)
    graphe=1;
end
if nargin<=2||isempty(nmax)
    nmax=1e7;
end
if nargout>=4
    untot=cell(1,Ntot);
end
un=zeros(1,Ntot);
p=un;
for i=1:Ntot
    u=i;
    test=1;
    k=0;
    if nargout>=4
        res=u;
    end
    while(test)
        k=k+1;
        if mod(u,2)==0
            u=u/2;
        else
            u=3*u+1;
        end
        if nargout>=4
            res=[res,u];
        end
        test=(u~=1)&&(k<nmax);
    end
    p(i)=k;
    un(i)=u;
    if nargout>=4
        untot{i}=res;
    end
    if k==nmax
       disp('Attention, nmax atteint');
    end
    if u~=1
       disp('Attention, valeur finale différente de 1');
    end
end
if graphe
    clf;
    plot(1:Ntot,p);
end
pmax=max(p);