Tech things and such
Article / Note
2017/04/22

Structures de Turing et modèle Gray-Scott de réaction-diffusion

Simulation interactive d'un célèbre système dynamique, tel que l'avait notamment imaginé Alan Turing dans The Chemical Basis of Morphogenesis pour expliquer comment certaines formes ou motifs complexes peuvent apparaître (par exemple, le pelage tacheté d'un léopard). Il s'agit ici d'un système à réaction-diffusion de Gray-Scott.






Clic et déplacer pour perturber le champ.

Quelques explications

On simule ici le comportement de deux substances chimiques. La première (A) est sécrétée dans le milieu à un certain débit (paramètre "feed"). Cette substance diffuse dans le milieu à une certaine vitesse et stimule la sécrétion d'une seconde substance (B). Celle-ci a la propriété de "capturer" la substance A. B diffuse dans le milieu à une plus grande vitesse que A. B est dégradé avec le temps et disparaît selon un taux (paramètre kill)

Ceci est modélisé par le système d'équations différentielles suivant :

$$ \frac{d}{dt}\left[ \begin{array}{c} A(x,y) \\\\ B(x,y) \end{array} \right]=\; \left[ \begin{array}{c} Da\Delta ^{2}A\; -\; A B^{2}\; +f\left( 1-A \right) \\\\ Db\Delta ^{2}B\; +\; AB^{2}\; -\left( f+k \right)B \end{array} \right] $$

Avec :

$$ A,B = concentration\ de\ A,B\ au\ point\ (x,y)$$ $$ Da,Db = vitesse\ de\ diffusion\ de\ A,B $$ $$ \Delta ^{2}A,\Delta ^{2}B = Laplacien\ de\ A,B\ au\ point\ (x,y) $$ $$ f = taux\ de\ sécrétion\ de\ A $$ $$ k = taux\ de\ suppression\ de\ B $$

Toutes les combinaisons de paramètres ne fonctionnent pas nécessairement. Le système exhibe une très forte sensibilité aux conditions initiales. Quoiqu'il en soit, le rapport Da/Db devra être proche de 2 avec Da ≤ 0,5. k peut varier entre 0.01 et 0.1. Le paramètre k se situera entre 0,04 et 0,07.

La page de Robert Munafo propose un diagramme d'exploration de ces paramètres très intéressant.

Ces paramètres f et k pourraient également varier dans le plan, la forme du motif serait alors progressivement modifiée... Je mets ceci dans ma todo list !

Algorithme

La fonction ci-après est exécutée 20 fois par frame update. u représente la concentration de a, v celle de b. On dessine ensuite la valeur de u pour chaque pixel du canvas.


function turing() {
    for(var y = 0 ; y < H ; y++) {
        for(var x = 0 ; x < W ; x++) {
            var index = x + y * W;
            // concentrations contient les concentrations de u et v + leurs dérivées. 
            var values = concentrations[index];

            // calcule des voisins (avec sécurité pour les bords du canvas).
            var up      = y > 0 ?       concentrations[index-W] : {u:values.u,v:values.v};
            var down    = y < (H-1) ?   concentrations[index+W] : {u:values.u,v:values.v};
            var left    = x > 0 ?       concentrations[index-1] : {u:values.u,v:values.v};
            var right   = x < (W-1) ?   concentrations[index+1] : {u:values.u,v:values.v};

            var upLeft      = y > 0 && x > 0 ?          concentrations[index-W-1] : {u:values.u,v:values.v};
            var downLeft    = y < (H-1) && x > 0 ?      concentrations[index+W-1] : {u:values.u,v:values.v};
            var upRight     = y > 0 && x < (W-1) ?      concentrations[index+1-W] : {u:values.u,v:values.v};
            var downRight   = y < (H-1) && x < (W-1) ?  concentrations[index+1+W] : {u:values.u,v:values.v};

            // calcul des laplaciens nécessaire pour la diffusion. Très similaire avec un algo de edge sharp        
            var laplu = 0.2*up.u + 0.2*down.u + 0.2*right.u + 0.2*left.u + 0.05*upLeft.u + 0.05*upRight.u + 0.05*downLeft.u + 0.05*downRight.u - values.u;
            var laplv = 0.2*up.v + 0.2*down.v + 0.2*right.v + 0.2*left.v + 0.05*upLeft.v + 0.05*upRight.v + 0.05*downLeft.v + 0.05*downRight.v - values.v;

            // Calcul des taux de variations (dérivées) de u et v. On retrouve ici les deux équations
            // Cu et Cv sont les taux de diffusion de A et B.           
            values.du = Cu * laplu - values.u*values.v*values.v + f*(1.0-values.u);
            values.dv = Cv * laplv + values.u*values.v*values.v - (k+f)*values.v;

            // Update des concentrations
            values.u = values.u + values.du;
            values.v = values.v + values.dv;
        }
    }
}

Quelques références :

Ce travail est en grande partie basé sur les études approfondies du sujet exposées sur le site de Robert Munafo. Qu'il soit ici vivement remercié ! Voici un lien vers son excellent simulateur. La liste de paramètres prédéfinis "feed" et "kill" est un emprunt empreint (!) d'humilité et de reconnaissance.

>> Réagir à cet article