Comment installer et utiliser ADiMat pour Matlab ?

Qu'est-ce que ADiMat ?

ADiMat signifie Automatisches Differenzieren für Matlab, c'est donc un programme tournant sous Matlab qui permet de "différencier" (calculer la différentielle, ou "dériver") n'importe quelle fonction écrite en Matlab, autrement dit une bête boite noire :

function sorties = f(entrees)
    % le corps de la fonction
    % transformant les entrées en sorties
end

Pourquoi utiliser AdiMat ?

Parce que parfois la fonction peut être extrêmement compliquée, et donc la jacobienne très difficile à écrire voir impossible à expliciter. Nous verrons plus loin, après l'installation, un exemple illustrant un cas difficile et la facilité avec laquelle ADiMat fournit une solution.

Installation de ADiMat

1) Télécharger l'archive correspondant à votre système d'exploitation et l'extraire dans votre HOME (~/).

tar xzvf adimat-0.5.7-3269-GNU_Linux-x86_64.tar.gz

2)  Dans votre dossier d'installation Matlab (chez moi ~/Matlab/) ouvrez le fichier startup.m et ajoutez ces deux lignes à la fin :

setenv('ADIMAT_HOME', '~/adimat-0.5.7-3269-GNU_Linux-x86_64');
ADiMat_startup

Ce qui renseignera le chemin vers le dossier d'installation d'ADiMat via la variable d'environnement ADIMAT_HOME. Et la seconde ligne initialisera le programme ADiMat au démarrage de Matlab.

Parfois le fichier startup.m ne se trouve pas à la racine du dossier Matlab mais plus profondément dans l'arborescence, pour le localiser placez vous dans le dossier Matlab et tapez simplement :

locate startup.m

3) Placez le fichier de démarrage suivant se trouvant dans votre dossier d'installation ADiMat :

~/adimat-0.5.7-3269-GNU_Linux-x86_64/share/adimat/runtime/ADiMat_startup.m

dans le même dossier Matlab où se trouve startup.m

4) Relancez Matlab, et l'installation est terminée.

Utilisation d'ADiMat

C'est très simple, une fois définie votre fonction :

function sorties = f(entrees)
    % le corps de la fonction
    % transformant les entrées en sorties
end

Le calcul de jacobienne évaluée au point "(arg1,...,argN)" s'obtient alors tout simplement (cf. ici ) :

J = admDiffFor(@f, 1, arg1, ..., argN)

Pour plus de détails consultez la documentation

Matlab vous demande alors si vous l'autorisez à communiquer avec le serveur d'ADiMat, répondez "yes". Si tout se passe bien il devrait vous calculer la matrice jacobienne J. Il se peut que vous obteniez à la place des erreurs, par exemple pour ma part je n'avais pas l'interpréteur /bin/zsh d'installé, ni la librairie libssl0.9.8

sudo apt-get install zsh
sudo apt-get install libssl0.9.8

Exemple d'utilisation

Un exemple est déjà donné sur le site officiel, mais il est relativement simpliste et on pourrait calculer la jacobienne à la main, ce qui ne met pas en valeur la puissance d'ADiMat. Pour l'illustrer son intérêt je vais traiter un exemple qui m'a sauvé la mise dans le cadre d'un filtrage de Kalman.

On part d'une équation différentielle d'ordre 2 que l'on réduit à l'ordre 1 en considérant 2 variables regroupées dans le vecteur , ainsi l'équation différentielle se réécrit

Nous discrétisons ensuite l'équation différentielle suivant un schéma de Runge-Kutta d'ordre 4 (RK4) :

avec

Nous avons relié l'état à l'état suivant via une certaine fonction :

Dont l'implémentation est bien sûr :

function [z1 v1] = RK4(z,v)
global h cstar
k1 = model(1,[z,v],cstar);
k2 = model(1,[z+h/2*k1(1),v+h/2*k1(2)],cstar);
k3 = model(1,[z+h/2*k2(1),v+h/2*k2(2)],cstar);
k4 = model(1,[z+h*k3(1),v+h*k3(2)],cstar);

z1 = z + h/6*[k1(1)+2*k2(1)+2*k3(1)+k4(1)];
v1 = v + h/6*[k1(2)+2*k2(2)+2*k3(2)+k4(2)];
end

Calculer la jacobienne de la fonction à la main est difficile car comme vous pouvez le voir les sont des composées de la fonction avec elle-même. Ce qui signifie que pour calculer :

vous allez devoir calculer ces dérivées partielles en utilisant plusieurs fois (à 4 niveaux car est une composée de 4 fois) la règle de la chaîne, ce qui est, vous en conviendrez, assez ardu.

ADiMat va faire le travail à votre place en utilisant une technique s'appuyant également sur la règle de la chaîne, à savoir la différentiation automatique. Et en une ligne vous obtenez le résultat escompté :

JRK4(n) = admDiffFor(@RK4, 1, z(n), v(n));