diff --git a/Optimisation/TP/TP-diff-finies/C14_diff_finies.m b/Optimisation/TP/TP-diff-finies/C14_diff_finies.m
new file mode 100644
index 0000000000000000000000000000000000000000..aee8962b08b048456bff5abc981869b4fa169989
--- /dev/null
+++ b/Optimisation/TP/TP-diff-finies/C14_diff_finies.m
@@ -0,0 +1,210 @@
+%*********************************************************************************************
+% Fichier ~gergaud/ENS//Users/gergaud/ENS/Optim1a/TP-optim/TP-diff-finies/C14_diff_finies.m *
+% Décembre 2017 *
+% Université de Toulouse, INP-ENSEEIHT *
+%*********************************************************************************************
+%
+%
+% ======================================================================
+% Ce programme estime les paramètres dans le problème d'estimation
+% des paramètres pour le problème de la datation par le carbonne 14
+% Programme MATLAB
+% Le modèle est:
+% A(t)=A0exp(-lambda t)
+% On cherche beta = (A0,lambda) de façon à ce que le modèle
+% approche le mieux les données au sens des moindres carrés
+% Les programmes et algorithmes utilisés pour la minimisation sont:
+% Gauss-Newton et leastsq
+% =======================================================================
+%
+clear all; close all; format shorte
+% Initialisation
+%
+diary res_C14_diff_finies.txt
+% Données:
+
+
+Ti=[500 1000 2000 3000 4000 5000 6300];
+Ai=[14.5 13.5 12.0 10.8 9.9 8.9 8.0];
+% Visualisation dans IR^3 de la fonction f à minimiser
+xmin=9;
+xmax=20;
+nx=100;
+pasx=(xmax-xmin)/nx;
+ymin=0;
+% ymin=-0.0001;
+ymax=0.0003;
+ny=100;
+pasy=(ymax-ymin)/ny;
+
+x=xmin:pasx:xmax;
+y=ymin:pasy:ymax;
+[X,Y]=meshgrid(x,y);
+[m,n]=size(X);
+for i=1:m,
+ for j=1:n,
+ Z(i,j)=fun_ref(@(beta) res_C14_ref(beta,Ti,Ai),[X(i,j),Y(i,j)]);
+ end;
+end;
+figure(1)
+mesh(X,Y,Z)
+xlabel('A_0');
+ylabel('\lambda')
+zlabel('f(A_0,\lambda')
+
+% Visualisation des données et des courbes de niveaux
+
+figure(2)
+subplot(1,2,1) % Données
+hold on;
+plot(Ti,Ai,'o')
+rect=[0 7000 0 18]; % dimensions des axes
+axis(rect)
+subplot(1,2,2) % Courbes de niveaux
+hold on
+contour(X,Y,Z,50);
+xlabel('A_0')
+ylabel('\lambda')
+hold on
+%
+% Gauss-Newton
+% ------------
+% Point de départ
+beta0 = [10 0.0001]'; % Newton, fminunc et leastsq converge
+% beta0=[15 0.001]' % Newton, fminunc et leastsq divergent
+%beta0=[15 0.0005]' % Newton diverge, fminu et leastsq convergent
+%beta0=[10 0.0005]'
+True_Jacobienne = false; % false : approximation par différence finies
+% Liste pour le choix de la méthode de différences
+list_diff_finies = {'avants','centrees'};
+%%%%% MODIFIER LE NUMERO APRES AVOIR CODE LA FONCTION ASSOCIEE DANS LE FICHIER ASSOCIE %%%%%%%%%%%
+methode_finite_diff = list_diff_finies{1};
+
+
+% On trace la courbe de départ
+T = linspace(0,6300,100);
+A = beta0(1)*exp(-beta0(2)*T);
+figure(2)
+subplot(1,2,1)
+plot(T,A)
+figure(2)
+subplot(1,2,2)
+text(beta0(1),beta0(2),'o point de depart')
+% pause
+%
+disp('Algorithme de Gauss-Newton')
+if ~True_Jacobienne
+ ndigits = 9; % précision ndigits
+ fprintf("ndigits %d \n", ndigits)
+end
+disp( strcat('Vraie Jacobienne = ' , sprintf(' %d', True_Jacobienne) ) )
+disp( strcat('Methode diff finies = ' + string(methode_finite_diff)) ) ;
+
+disp('------------------------------------------------------------------')
+disp(' k ||f''(beta)|| f(beta) ||s|| exitflag ')
+disp('------------------------------------------------------------------')
+
+if True_Jacobienne
+ Jres = J_res_C14_ref(beta0, Ti);
+else
+ switch methode_finite_diff
+ case 'avants'
+ Jres = diff_finies_avant( @(beta) res_C14_ref(beta,Ti,Ai), beta0, ndigits);
+
+ case 'centrees'
+ Jres = diff_finies_centree( @(beta) res_C14_ref(beta,Ti,Ai), beta0, ndigits);
+ end
+end
+
+res = res_C14_ref(beta0,Ti,Ai);
+gradf_beta = norm(Jres'*res);
+f_beta = fun_ref(@(beta) res_C14_ref(beta,Ti,Ai), beta0);
+disp([0 gradf_beta f_beta]);
+Beta = [beta0];
+
+for i=1:10
+ option = [sqrt(eps) 1.e-12 i];
+ if True_Jacobienne
+ [beta, norm_gradf_beta, f_beta, res, norms, k, exitflag] = GN_ref(@(beta) res_C14_ref(beta,Ti,Ai), ...
+ @(beta) J_res_C14_ref( beta, Ti), beta0, option);
+ else
+ switch methode_finite_diff
+ case 'avants'
+ [beta, norm_gradf_beta, f_beta, res, norms, k, exitflag] = GN_ref(@(beta) res_C14_ref(beta,Ti,Ai), ...
+ @(beta) diff_finies_avant(@(beta) res_C14_ref(beta,Ti,Ai), beta, ndigits), beta0, option);
+ case 'centrees'
+ [beta, norm_gradf_beta, f_beta, res, norms, k, exitflag] = GN_ref(@(beta) res_C14_ref(beta,Ti,Ai), ...
+ @(beta) diff_finies_centree(@(beta) res_C14_ref(beta,Ti,Ai), beta, ndigits), beta0, option);
+ end
+ end
+ Beta = [Beta beta];
+ disp([k norm_gradf_beta f_beta norms exitflag])
+% Visualisation
+ A=beta(1)*exp(-beta(2)*T);
+ subplot(1,2,1)
+ plot(T,A)
+% eval(['print -depsc fig_GN_courbe' int2str(i) '_C14'])
+ subplot(1,2,2)
+ plot(beta(1),beta(2),'o')
+end
+disp('------------------------------------------------------------------')
+disp('betak')
+disp(Beta);
+figure(1)
+
+print -depsc C14_fig_Gauss_Newton1
+figure(2)
+
+print -depsc C14_fig_Gauss_Newton2
+diary
+
+
+% ======================================================================
+% Calcul de la fonction à minimiser
+% =======================================================================
+function sce=fun_ref(res,beta)
+% Input:
+% ------
+% res : fonction résidus
+% res : \IR^p --> \IR^n
+% beta : vecteur des paramètres
+% real(p)
+%
+% Output:
+% -------
+% sce : (1/2)||res(beta)||^2
+% real
+%
+% Calcul des résidus
+res = res(beta);
+res = res(:);
+sce=0.5*res'*res;
+end
+%
+% ======================================================================
+% "résidus" pour le probleme de datation par le carbonne 14
+% =======================================================================
+%
+%
+function r=res_C14_ref(beta,Ti,Ai)
+% Calcul des résidus
+%
+r=Ai-beta(1)*exp(-beta(2)*Ti);
+r=r(:);
+end
+%
+% ======================================================================
+% Cette fonction calcule la derivee de la fonction residu pour l'exemple
+% datation par le carbonne 14
+% =======================================================================
+%
+%
+function J_res=J_res_C14_ref(beta,Ti)
+%
+%
+Ti = Ti(:);
+J_res = [-exp(-beta(2)*Ti) (beta(1)*exp(-beta(2)*Ti)).*Ti];
+end
+
+
+
diff --git a/Optimisation/TP/TP-diff-finies/GN_ref.p b/Optimisation/TP/TP-diff-finies/GN_ref.p
new file mode 100644
index 0000000000000000000000000000000000000000..db427cc05ccf407569d1ecf981b8640159b04529
Binary files /dev/null and b/Optimisation/TP/TP-diff-finies/GN_ref.p differ
diff --git a/Optimisation/TP/TP-diff-finies/TP_diff_finies.pdf b/Optimisation/TP/TP-diff-finies/TP_diff_finies.pdf
new file mode 100644
index 0000000000000000000000000000000000000000..5f798a2fc2d9ea6c72ed4df9a9c7df14d4613d76
Binary files /dev/null and b/Optimisation/TP/TP-diff-finies/TP_diff_finies.pdf differ
diff --git a/Optimisation/TP/TP-diff-finies/diff_finies_avant.m b/Optimisation/TP/TP-diff-finies/diff_finies_avant.m
new file mode 100644
index 0000000000000000000000000000000000000000..5f1491ada3ae3ac999da6bd9830e66d711e73301
--- /dev/null
+++ b/Optimisation/TP/TP-diff-finies/diff_finies_avant.m
@@ -0,0 +1,38 @@
+% Auteur : J. Gergaud
+% décembre 2017
+% -----------------------------
+%
+
+
+
+function Jac = diff_finies_avant(fun,x,option)
+%
+% Cette fonction calcule les différences finies avant sur un schéma
+% Paramètres en entrées
+% fun : fonction dont on cherche à calculer la matrice jacobienne
+% fonction de IR^n à valeurs dans IR^m
+% x : point où l'on veut calculer la matrice jacobienne
+% option : précision du calcul de fun (ndigits)
+%
+% Paramètre en sortie
+% Jac : Matrice jacobienne approximé par les différences finies
+% real(m,n)
+% ------------------------------------
+Jac = 0;
+end
+
+function s = sgn(x)
+% fonction signe qui renvoie 1 si x = 0
+if x==0
+ s = 1;
+else
+ s = sign(x);
+end
+end
+
+
+
+
+
+
+
diff --git a/Optimisation/TP/TP-diff-finies/diff_finies_centree.m b/Optimisation/TP/TP-diff-finies/diff_finies_centree.m
new file mode 100644
index 0000000000000000000000000000000000000000..01dc995fcdaa9d6159483eb489dcab00a92c9eb0
--- /dev/null
+++ b/Optimisation/TP/TP-diff-finies/diff_finies_centree.m
@@ -0,0 +1,36 @@
+% Auteur : J. Gergaud
+% décembre 2017
+% -----------------------------
+%
+
+
+
+function Jac= diff_finies_centree(fun, x, option)
+%
+% Cette fonction calcule les différences finies centrées sur un schéma
+% Paramètres en entrées
+% fun : fonction dont on cherche à calculer la matrice jacobienne
+% fonction de IR^n à valeurs dans IR^m
+% x : point où l'on veut calculer la matrice jacobienne
+% option : précision du calcul de fun (ndigits)
+%
+% Paramètre en sortie
+% Jac : Matrice jacobienne approximé par les différences finies
+% real(m,n)
+% ------------------------------------
+Jac = 0;
+end
+
+function s = sgn(x)
+% fonction signe qui renvoie 1 si x = 0
+if x==0
+ s = 1;
+else
+ s = sign(x);
+end
+end
+
+
+
+
+
diff --git a/Optimisation/TP/TP-diff-finies/exem1_diff_finies.m b/Optimisation/TP/TP-diff-finies/exem1_diff_finies.m
new file mode 100644
index 0000000000000000000000000000000000000000..6fcea57406a0cbd5d25c5405d2d6537c467fd8d9
--- /dev/null
+++ b/Optimisation/TP/TP-diff-finies/exem1_diff_finies.m
@@ -0,0 +1,173 @@
+%****************************************************************************************************
+% Fichier ~gergaud/ENS//Users/gergaud/ENS/Optim1a/TP-optim/TP-diff-finies/exem1_diff_finies_ref.m *
+% Décembre 2017 *
+% Université de Toulouse, INP-ENSEEIHT *
+%****************************************************************************************************
+%
+%
+% ======================================================================
+% Ce programma résout par l'algorithme de Gauss-Newton le problème
+% d'optimisation
+% Min f(x) = (1/2)((x_1^2-x_2)^2+(1-x_1)^2)=(1/2)||r(x)||^2
+% =======================================================================
+%
+
+
+
+clear all; close all; format shorte
+diary res_exem1_diff_finies.txt
+
+
+
+
+
+% Initialisation
+%
+% Visualisation dans IR^3 de la fonction f à minimiser
+xmin=-1;
+xmax=2;
+nx=100;
+pasx=(xmax-xmin)/nx;
+ymin=-1;
+% ymin=-0.0001;
+ymax=4;
+ny=100;
+pasy=(ymax-ymin)/ny;
+x=xmin:pasx:xmax;
+y=ymin:pasy:ymax;
+[X,Y]=meshgrid(x,y);
+[m,n]=size(X);
+for i=1:m
+ for j=1:n
+ Z(i,j)=fun_ref([X(i,j),Y(i,j)]);
+ end
+end
+figure(1)
+mesh(X,Y,Z)
+xlabel('x_1');
+ylabel('x_2')
+zlabel('f(x)')
+
+% Visualisation des données et des courbes de niveaux
+
+figure(2) % Courbes de niveaux
+
+hold on
+contour(X,Y,Z,50);
+xlabel('x_1')
+ylabel('x_2')
+hold on
+%
+% Gauss-Newton
+% ------------
+% Point de départ
+beta0 = [500 100]';
+text(beta0(1),beta0(2),'o point de depart')
+True_Jacobienne = false; % false : approximation par différence finies
+% Liste pour le choix de la méthode de différences
+list_diff_finies = {'avants','centrees'};
+%%%%% MODIFIER LE NUMERO APRES AVOIR CODE LA FONCTION ASSOCIEE DANS LE FICHIER ASSOCIE %%%%%%%%%%%
+methode_finite_diff = list_diff_finies{1};
+
+
+
+disp('Algorithme de Gauss-Newton')
+disp( strcat('Vraie Jacobienne =' , sprintf('%d', True_Jacobienne) ) );
+disp( strcat('Methode diff finies = ' + string(methode_finite_diff)) ) ;
+
+if ~True_Jacobienne
+ ndigits = 8; % précision ndigits
+ fprintf("ndigits : %d \n", ndigits)
+end
+
+disp('------------------------------------------------------------------')
+disp(' k ||f''(beta)|| f(beta) ||delta|| exitflag ')
+disp('------------------------------------------------------------------')
+
+if True_Jacobienne
+ Jres = J_res1(beta0);
+else
+ switch methode_finite_diff
+ case 'avants'
+ Jres = diff_finies_avant(@res1,beta0, ndigits);
+ case 'centrees'
+ Jres = diff_finies_centree(@res1,beta0, ndigits);
+ end
+end
+
+res = res1(beta0);
+gradf_beta = norm(Jres'*res);
+f_beta = fun_ref(beta0);
+disp([0 gradf_beta f_beta]);
+Beta = [beta0];
+for i=1:2
+ option = [sqrt(eps) 1.e-12 i];
+ if True_Jacobienne
+ [beta, norm_gradf_beta, f_beta, res, norms, k, exitflag] = GN_ref(@res1, @J_res1, beta0, option);
+ else
+
+ switch methode_finite_diff
+ case 'avants'
+ [beta, norm_gradf_beta, f_beta, res, norms, k, exitflag] = GN_ref(@res1, @(beta) diff_finies_avant(@res1,beta,ndigits), beta0, option);
+ case 'centrees'
+ [beta, norm_gradf_beta, f_beta, res, norms, k, exitflag] = GN_ref(@res1, @(beta) diff_finies_centree(@res1,beta,ndigits), beta0, option);
+ end
+end
+ Beta = [Beta beta];
+ disp([k norm_gradf_beta f_beta norms exitflag])
+% Visualisation
+ plot(beta(1),beta(2),'o')
+end
+disp('------------------------------------------------------------------')
+disp('betak')
+disp(Beta);
+figure(1)
+print -depsc exem_1_fig_GN1_exem1
+figure(2)
+print -depsc exem_1fig_GN2_exem1
+diary
+
+
+% ======================================================================
+% Calcul de la fonction à minimiser
+% =======================================================================
+
+
+
+% ======================================================================
+% "résidus" pour le probleme:
+% Min f(x)=0.5*((x_1^2-x_2)^2+(1-x_1)^2)
+% =======================================================================
+%
+%
+function r=res1(x)
+r=[x(1)^2-x(2)
+ 1-x(1)];
+end
+
+
+% ======================================================================
+% Calcul de la fonction objective à minimiser
+% Min f(x)=0.5*((x_1^2-x_2)^2+(1-x_1)^2)
+% =======================================================================
+%
+%
+function sce=fun_ref(x)
+r=[x(1)^2-x(2)
+ 1-x(1)];
+ sce = 0.5*r'*r;
+end
+
+% ========================================================================
+% Cette fonction calcule la derivee de la fonction residu pour l'exemple 1
+% ========================================================================
+%
+%
+function J_res=J_res1(x)
+%
+%
+J_res = [2*x(1) -1 ; -1 0];
+end
+
+
+
diff --git a/Optimisation/TP/TP-diff-finies/main_diff_finies_cosinus.m b/Optimisation/TP/TP-diff-finies/main_diff_finies_cosinus.m
new file mode 100644
index 0000000000000000000000000000000000000000..294112c47ad15fb150a2d2888ebfeabfdef91ac2
--- /dev/null
+++ b/Optimisation/TP/TP-diff-finies/main_diff_finies_cosinus.m
@@ -0,0 +1,211 @@
+% Auteur : Joseph Gergaud, INP-ENSEEIHT & IRIT
+% Date : décembre 2017
+%
+% Test de l'algorithme de différences finies avant sur la fonction cosinus
+% fun1(x) = cos(x)
+% fun2(x) = cos(x) + alpha*rand(1)
+% Avec 2 points différents
+%
+
+% On réinitialise l'environnement
+%
+clear;
+close all;
+clc;
+
+% Des paramètres pour l'affichage
+%
+tirets = '------------------------------------------';
+LW = 1.5;
+set(0,'defaultaxesfontsize' , 14 , ...
+ 'DefaultTextVerticalAlignment' , 'bottom', ...
+ 'DefaultTextHorizontalAlignment', 'left' , ...
+ 'DefaultTextFontSize' , 12 , ...
+ 'DefaultFigureWindowStyle' ,'docked');
+
+% On choisit un format long pour plus de détails sur les résultats numériques
+%
+format shortE;
+
+% Exercice 2.2
+%-------------------------------------------------------------------------------------------------------------
+% Test de l'algorithme de différences finies avant sur la fonction cosinus
+%
+disp(tirets);
+disp('Test de l''algorithme de différences finies sur la fonction cosinus');
+
+% Les differentes fonctions et la jacobienne theorique
+fun1 = @(x) -cos(x);
+fun2 = @(x) -cos(x) + 1.e-8*rand(1);
+true_jac = @(x) sin(x);
+
+% Variables pour le choix de la methode de difference finie
+list_diff_finies = {'avants','centrees','complexes'};
+%%%%% MODIFIER LE NUMERO APRES AVOIR CODE LA FONCTION ASSOCIEE %%%%%%%%%%%
+methode_finite_diff = list_diff_finies{1};
+
+% Calcul de l'erreur sur la fonction dans le cas bruité et non bruité
+omega_true_fun = eps((1));
+omega_per_fun = 1e-8;
+
+
+
+% Estimation du h* en fonction de la difference finie choisie
+switch methode_finite_diff
+ case 'avants'
+ h_star_true_fun = sqrt(omega_true_fun);
+ h_star_per_fun = sqrt(omega_per_fun);
+ case 'centrees'
+ h_star_true_fun = omega_true_fun^(1/3.);
+ h_star_per_fun = omega_per_fun^(1/3.);
+ otherwise
+ h_star_true_fun = sqrt(omega_true_fun);
+ h_star_per_fun = sqrt(omega_true_fun);
+end
+
+
+
+
+
+% Points ou l'on souhaite effectuer les tests
+x0 = pi/3;
+x1 = 1.e6*pi/3;
+
+% Ordres pour faire les tests (16 + celui de h*)
+ordrestrue = ([1:floor(-log10(h_star_true_fun)), -log10(h_star_true_fun), ceil(-log10(h_star_true_fun)):16]);
+ordresper = ([1:floor(-log10(h_star_per_fun)), -log10(h_star_per_fun), ceil(-log10(h_star_per_fun)):16]);
+% Initialisation des vecteurs d'erreur
+
+err_x0 = zeros(1,length(ordrestrue));
+err_x0p = zeros(1,length(ordresper));
+err_x1 = zeros(1,length(ordrestrue));
+
+% Définition du schéma à employer
+
+
+switch methode_finite_diff
+ % Differences finies avant
+ case 'avants'
+ diff_finies = @(fun,x,h) forwardfiniteDiff(fun,x,h);
+ % Differences finies avant
+ case 'centrees'
+ diff_finies = @(fun,x,h) centredfiniteDiff(fun,x,h);
+ % Differences finies avant
+ case 'complexes'
+ diff_finies = @(fun,x,h) derivee_complexe(fun,x,h);
+end
+
+for i = 1:length(ordrestrue)
+ h = 10^(-ordrestrue(i));
+ err_x0(i) = abs(diff_finies(fun1,x0,h) - (true_jac(x0)));
+ err_x1(i) = abs(diff_finies(fun1,x1,h) - (true_jac(x1)));
+
+ h = 10^(-ordresper(i));
+ err_x0p(i) = abs(diff_finies(fun2,x0,h) - (true_jac(x0)));
+
+
+end
+
+% Affichage de la courbe pour x0
+affichage_erreur(ordrestrue, err_x0, h_star_true_fun, 'x0', methode_finite_diff, LW);
+
+% Affichage de la courbe pour x1
+affichage_erreur(ordrestrue, err_x1, h_star_true_fun, 'x1', methode_finite_diff, LW);
+% Affichage de la courbe pour x0p
+affichage_erreur(ordresper, err_x0p, h_star_per_fun, 'x0p', methode_finite_diff, LW);
+
+
+
+
+
+
+% Fonctions ---------------------------------------------------------------
+
+function affichage_erreur(ordres,err,h_star,x_str,methode_str,LW)
+ figure;
+ non_nul_element = find(err);
+ err = err(non_nul_element);
+ ordres = ordres(non_nul_element);
+ % Courbe des erreurs pour les differents ordres en bleu
+ loglog(10.^(-ordres), err, 'b', 'LineWidth', LW);
+ hold on;
+ % Ligne verticale pour situer le pas optimale h* en rouge
+ line([h_star, h_star], [min(err), max(err)], 'Color', 'r');
+ grid on;
+ xlabel('h');
+ ylabel('erreurs');
+
+ titre = ['Erreur des differences finies ' methode_str];
+
+
+
+ switch x_str
+ case 'x0'
+ titre = [titre, sprintf('\n sur la fonction cosinus en x_0= %s/3', '\pi')];
+ case 'x1'
+ titre = [titre, sprintf('\n sur la fonction cosinus en x_1= %s/3', '1.e6*\pi')];
+ case 'x0p'
+ titre = [titre, sprintf('\n sur la fonction cosinus perturbee en x_0= %s/3', '\pi')];
+ end
+
+ title(titre, 'HorizontalAlignment', 'center');
+ legend({'Numerique','Optimal'},'Location','SouthWest');
+ grid on;
+ drawnow;
+ % Enregistrement de l'image
+ print(['fig_diff_finies_' x_str],'-dpng');
+end
+
+
+function Jac = forwardfiniteDiff(fun,x,h)
+%
+% Cette fonction calcule les différences finies
+% Paramètres en entrées
+% fun : fonction dont on cherche à calculer la matrice jacobienne
+% fonction de IR^n à valeurs dans IR^m
+% x : point où l'on veut calculer la matrice jacobienne
+% h : pas
+%
+% Paramètre en sortie
+% Jac : Matrice jacobienne approximé par les différences finies
+% real(m,n)
+% ------------------------------------
+Jac = 0;
+end
+
+function Jac = centredfiniteDiff(fun, x, h)
+%
+%Cette fonction calcule les différences finies à l'aide d'un schéma centré
+%Paramètres en entrées
+% fun : fonction dont on cherche à calculer la matrice jacobienne
+% fonction de IR^n à valeurs dans IR^m
+% x : point où l'on veut calculer la matrice jacobienne
+% h : pas
+%
+% Paramètre en sortie
+% Jac : Matrice jacobienne approximé par les différences finies
+% real(m,n)
+% ------------------------------------
+Jac = 0;
+end
+
+function Jac = derivee_complexe(fun, x,h)
+
+%
+%Cette fonction calcule les différences finies à l'aide d'un schéma complexe
+%Paramètres en entrées
+% fun : fonction dont on cherche à calculer la matrice jacobienne
+% fonction de IR^n à valeurs dans IR^m
+% x : point où l'on veut calculer la matrice jacobienne
+% h : pas
+%
+% Paramètre en sortie
+% Jac : Matrice jacobienne approximé par les différences finies
+% real(m,n)
+% ------------------------------------
+Jac = 0;
+end
+
+
+
+