3. Phase réactionnelle

On considère ici, comme d’autres auteurs [2, 3, 4], que les groupements chlorés du TMC sont hydrolysés.

_images/apa_reactants.png

Fig. 4 Les deux monomères initiaux, numérotés. L’ordre des atomes est celui donné par le code smile.

Quand la réaction amide+acide a lieu, le type de l’atome N1 passe de amine à amide. Idem pour le type de l’atome C1 de l’acide carboxylique, etc. Ces changements de type doivent être pris en compte lors de la simulation. Pour cela, on doit créer un patron pre-réactif et un patron post-réactif qui vont permettre de définir les types des atomes impliqués avant et après réaction. Ces patrons ne doivent contenir que les atomes impliqués dans des liaisons, angles, dièdres, etc. qui sont susceptibles de changer ainsi les atomes de bord (edge atoms), les premiers atomes non-impliqués. C’est le rôle des fichiers de patron (template) moléculaire de préciser tous les changements qui interviennent au cours de la réaction. Compte-tenu de la topologie des réactifs, on doit donc distinguer les types de réactions suivantes :

  • mpd (0 substitution) + tmc (0 substitution) (notée mpd0s_tmc0s)

  • mpd (0 substitution) + tmc (1 substitution)

  • mpd (0 substitution) + tmc (2 substitution)

  • mpd (1 substitution) + tmc (0 substitution)

  • mpd (1 substitution) + tmc (1 substitution)

  • mpd (1 substitution) + tmc (2 substitution) (mpd1s_tmc2s)

Cet ensemble de templates permet de décrire toutes les réactions possibles dans ce matériau.

3.1. Le patron mpd0s_tmc0s

On présente ci-dessous Fig. 5 un patron dont peuvent être dérivés tous les autres patrons réactifs. Il représente la réaction entre un monomère de MPD et un monomère de TMC.

_images/mpd0s_tmc0s.png

Fig. 5 Le patron des 2 réactifs pour la réaction de type 1, mpd+tmc. En jaune les atomes “réactifs”, qui vont former une nouvelle liaison covalente, en rouge, les atomes qui seront éliminés lors de la réaction, en bleu les atomes de bord. Les atomes grisés représentent les parties non prises en compte explicitement dans le template car elles contiennent des atomes qui sont séparés de plus de 4 liaisons des atomes réactifs. Ces atomes ne sont donc pas numérotés.

3.2. Les autres patrons

Pour les templates qui concernent les réactions avec le TMC 0, 1 ou 2 fois substitués (i.e. ayant réagi 0, 1 ou 2 fois), il faut remplacer le type des atomes 19 et 22 des fonctions acides par le type amide, soit C470 par C234. Pour les templates qui concernent le MPD 0 ou 1 fois substitué (i.e. ayant réagi 0 ou 1 fois), il faut remplacer le type des atomes 1, 2 et 3 de la fonction amine par les types amide, soit N900 par N265, H909 par H241 et C916 par C266.

La numérotation des atomes ne changeant pas selon la réaction considérée, on peut utiliser le même fichier de mapping, mpd_tmc.map.

3.3. Mise en œuvre avec lammps

La première étape est l’addition dans le fichier topology.dat de tous les templates pre- et post- réactionels nécessaires. Comme indiqué ci-dessus, il en existe 6 de chaque type.

Ces templates sont ensuite lus par lammps et associés à un nom de molécule. Comme de nouvelles réactions vont avoir lieu qui risquent potentiellement d’ajouter de nouvelles liaisons, angles, dièdres, etc., LAMMPS demande de prévoir le redimensionnement des tableaux associés aux listes qui gèrent ces interactions.  :

read_data data.in   extra/special/per/atom    30 &
                    extra/bond/per/atom       40 &
                    extra/angle/per/atom      40 &
                    extra/dihedral/per/atom   40 &
                    extra/improper/per/atom   40

molecule mpd mpd.tmpl offset 0 0 0 0 0
molecule tmc tmc.tmpl offset 0 0 0 0 0
molecule mpd0s_tmc0s_pre  mpd0s_tmc0s_pre.tmpl  offset 0 0 0 0 0
molecule mpd0s_tmc0s_post mpd0s_tmc0s_post.tmpl offset 0 0 0 0 0
molecule mpd0s_tmc1s_pre  mpd0s_tmc1s_pre.tmpl  offset 0 0 0 0 0
molecule mpd0s_tmc1s_post mpd0s_tmc1s_post.tmpl offset 0 0 0 0 0
molecule mpd0s_tmc2s_pre  mpd0s_tmc2s_pre.tmpl  offset 0 0 0 0 0
molecule mpd0s_tmc2s_post mpd0s_tmc2s_post.tmpl offset 0 0 0 0 0
molecule mpd1s_tmc0s_pre  mpd1s_tmc0s_pre.tmpl  offset 0 0 0 0 0
molecule mpd1s_tmc0s_post mpd1s_tmc0s_post.tmpl offset 0 0 0 0 0
molecule mpd1s_tmc1s_pre  mpd1s_tmc1s_pre.tmpl  offset 0 0 0 0 0
molecule mpd1s_tmc1s_post mpd1s_tmc1s_post.tmpl offset 0 0 0 0 0
molecule mpd1s_tmc2s_pre  mpd1s_tmc2s_pre.tmpl  offset 0 0 0 0 0
molecule mpd1s_tmc2s_post mpd1s_tmc2s_post.tmpl offset 0 0 0 0 0

On peut ensuite utiliser le fix bond/react pour réaliser toutes les réactions nécessaires :

fix rxn all bond/react stabilization yes nptrxn 0.01 reset_mol_ids yes &
  react rxn1 all 1 0.0 5.0 mpd0s_tmc0s_pre  mpd0s_tmc0s_post  mpd_tmc.map &
  react rxn2 all 1 0.0 5.0 mpd0s_tmc1s_pre  mpd0s_tmc1s_post  mpd_tmc.map &
  react rxn3 all 1 0.0 5.0 mpd0s_tmc2s_pre  mpd0s_tmc2s_post  mpd_tmc.map &
  react rxn4 all 1 0.0 5.0 mpd1s_tmc0s_pre  mpd1s_tmc0s_post  mpd_tmc.map &
  react rxn5 all 1 0.0 5.0 mpd1s_tmc1s_pre  mpd1s_tmc1s_post  mpd_tmc.map &
  react rxn6 all 1 0.0 5.0 mpd1s_tmc2s_pre  mpd1s_tmc2s_post  mpd_tmc.map

On peut suivre l’avancement de la réaction grace aux variables associées au fix bond/react. Plutôt que de sauver ces variables dans le fichier de données thermodynamiques, ont crée un fichier dedié à ces informations, le fichier de nom générique rxn.out et on supprime l’affichage à l’écran/log :

variable time   equal step
variable r1   equal f_rxn[1]
variable r2   equal f_rxn[2]
variable r3   equal f_rxn[3]
variable r4   equal f_rxn[4]
variable r5   equal f_rxn[5]
variable r6   equal f_rxn[6]
variable sum  equal f_rxn[1]+f_rxn[2]+f_rxn[3]+f_rxn[4]+f_rxn[5]+f_rxn[6]

fix rxnfile all print ${perthermo} " ${time} ${sum} ${r1} ${r2} ${r3} ${r4} ${r5} ${r6} " file rxn.out screen no

Enfin, on lance la simulation avec un fix npt agissant uniquement sur les atomes qui ne sont pas impliqués dans la réaction et en annulant le moment d’ensemble. La structure finale est sauvegardée dans le fichier de nom générique data.out et les informations de champ de force et du fix bond/react sont retirées :

fix zeromom nptrxn_REACT momentum ${permom} linear 1 1 1 angular
...
fix 1 nptrxn_REACT npt temp 298.0 298.0 100.0 iso 0.986923 0.986923 1000
...
write_data data.out nocoeff nofix

3.4. Tests et essais du ‘’fix bond/react’’

On présente ci-dessous les premiers résultats obtenus avec le fix bond/react sans imposer de conditions particulières sur les conditions de réaction, outre la distance maximale de 5 Å entre sites réactifs.

                       025        150        500
---------------------------------------------------
MPD/TMC             300/1200   300/200    1000/200
Mol kept               234          4        122
Ats/largest mol        145       7175        456
---------------------------------------------------
mpd_0 + tmc_0          300        195        200
mpd_0 + tmc_1            0         79        178
mpd_0 + tmc_2            0         23        134
mpd_1 + tmc_0          230          5          0
mpd_1 + tmc_1           68        111         22
mpd_1 + tmc_2            2        129         66
Total rxn              600        542

Tous ces systèmes peuvent former un maximum de 600 fonctions amides, mais le choix de la concentration initiale influe beaucoup sur la nature du réseau formé. Par exemple, il n’y a eu que 2 réactions entre un mpd une fois substitué (mpd_1) et un tmc 2 fois substitué (tmc_2), ce qui conduit à un réseau faiblement réticulé et une taille de la plus grosse molécule extrêmement faible. Dans le mélange MPD/TMC=1.5, on converge vers 540 fonctions amides mais 129 réactions mpd_1/tmc_2, donc un système beaucoup plus réticulé avec une molécule de 7175 atomes, puis de 62 puis des 2 dimères et 3 monomères de mpd et un total de 7 molécules.

max_rxn Afin de modifier la structure du réseau on peut imposer le nombre maximum de réactions par type.

Ordre des réactions dans le fix L’ordre des réactions dans le fix bond/react n’a aucune influence : sur les 6 réactions impliquées, on note une seule différence liée vraisemblablement à la sensibilité aux conditions initiales.