function [er,X]=sqrtm_iteration_newton(A,X0,epsilon)

% 
%  sqrtm_iteration_newton : calcul d'une racine carrée de matrice par la méthode de Newton.
%
%
% *********************************************************
%
%  [er,X]=sqrtm_iteration_newton(A,X0,espilon)
%      calcule les différents itérés de
%      la méthode de Newton pour la recherche d'une racine de matrice.
%      A chaque itération, elle affiche et stocke le résultat X_k ainsi que l'erreur
%         er(k)=|||X_k^2-A|||.
%	
%	variables d'entrées :
%    * A :matrice dont on veut une racine.
%    * X0 : matrice initiale
%    * epsilon : paramètre de précision
%	variables de sortie :
%    * er: suite des erreurs, définies par er(k)=norm(X_k^2-A).
%    * X : suite des itérés (matriciels) calculés.
% 
% ATTENTION : nombre maximum d'itération défini par nmax
% 
% ************ Fonctions auxiliaires utilisées ************
%
%       aucune
%
% *********************************************************
%

% nombre maximum d'itération défini par nmax
nmax=100;

% Contrôles d'entrée
[n,p]=size(A);
if (n~=p) 
   error('attention, A n''est pas carrée : arrêt du programme');
end
if (norm(A*X0-X0*A,inf)~=0)
   disp('attention, la matrice initiale et A ne commutent pas !!');
end

   

% Corps d'algorithme
er=[];
x=X0;
X=x;
k=0;
test=(norm(X0^2-A,inf)~=0);
while(test)
   y=(x+A/x)/2;
   disp(['itération numéro ',int2str(k)]);
   y2=y^2;
   erloc=norm(y2-A);
   er=[er erloc];
   X=[X y];
   test=~((norm(y2-A,inf)==0)|(k>=nmax)|(norm(x-y,inf)<=epsilon));
   disp('itéré : ');
   disp(y);
   disp(['erreur : ',num2str(erloc)]);
   disp('appuyez sur une touche pour continuer');
   pause;
   k=k+1;
   x=y;
end
if (n>=nmax)
   disp('attention, n>nmax');
end
