Modélisation moléculaire sur Raspberry Pi

Bonjour à tous,

Ouvrons une petite discussion sur la chimie computationnelle!

J'ai installé ce week-end GAMESS sur mon Raspberry Pi (modèle B). J'ai plus ou moins suivi cette procédure mais tout ne s'est pas bien déroulé… Au final, après pas mal de bidouillage, ça tourne. Malheureusement, comme prévu, c'est très lent: 25.6 minutes pour optimiser la géométrie de CH4 en M06-L/DZV… Moi qui voudrais faire des complexes métalliques… Voici mon input, pour ceux qui seraient intéressés:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16

 $contrl runtyp=optimize dfttyp=M06-L $end
 $basis gbasis=DZV $end

 $scf dirscf=.true. $end
 $statpt opttol=0.0005 nstep=50 $end
 $system mwords=50 $end

 $data
Methane
C1
C    6.0    -5.29    -2.22    0.00
H    1.0    -4.22    -2.22    0.00
H    1.0    -5.65    -2.71    0.88
H    1.0    -5.65    -2.71    -0.86
H    1.0    -5.65    -1.21    -0.01
 $end

(je me suis inspiré de ce blog pour les paramètres de statpt et system.

Malgré tout je pense que ça a du potentiel. Si on est patient, on peut faire des belles choses sans mobiliser son ordinateur et en consommant peu d'électricité.

D'autres retours sur GAMESS, Rapsberry Pi, … ?

Des retours sur Raspberry Pi, non, j'avoue. Des retours sur GAMESS, par contre, oui, un certainS nombreS (mais je pense que ça embêterai les gens plus qu'autre chose ). Dans tout les cas, si j'ai un jour le temps et l'énergie, j'écrirai un tuto de chimie quantique, et vu que GAMESS est gratuit (même si sous licence), je pense proposer de s'en servir pour faire deux ou trois trucs drôle. Mais … Euh … Pas tout de suite

Dans tout les cas, je trouve l'idée intéressante. Malgré la petite taille du machin et le processeur pas forcément efficace, ça fait toujours plaisir de ne pas manger un cœur de sa machine perso pour faire tourner GAMESS (ou plus). En plus, question refroidissement, tu peux faire des trucs un minimum intéressants, vu la surface à refroidir.

Moi, je dis "pourquoi pas". Et, bien entendu, "je suis jaloux"

J'arrive pas à me rendre compte de ce que représente 25 minutes pour ce calcul puisque je ne suis pas familier avec ces méthodes. Tu pourrais donner quelques points de comparaison ? Par exemple un ordinateur conventionnel, et une simulation sur un node de super-ordinateur ?

J'avais pensé faire ça pour VASP, mais j'avais peur que les limitations en mémoire vive ne permette pas de le faire tourner, indépendamment du temps qu'on lui laisse. Tu penses que ça marchera pour de plus grosses molècules ?

Dans tout les cas, si j'ai un jour le temps et l'énergie, j'écrirai un tuto de chimie quantique, et vu que GAMESS est gratuit (même si sous licence), je pense proposer de s'en servir pour faire deux ou trois trucs drôle.

Un autre code, open-source pour sa part : http://www.quantum-espresso.org/. Ça fonctionne plutôt pas mal, j'ai un petit projet de métadynamique ab initio avec et ça tourne bien.

J'arrive pas à me rendre compte de ce que représente 25 minutes pour ce calcul puisque je ne suis pas familier avec ces méthodes. Tu pourrais donner quelques points de comparaison ? Par exemple un ordinateur conventionnel, et une simulation sur un node de super-ordinateur ?

Pour la simulation en question, je prends 1s par calcul SCF sur 8 cœurs, pour 5 atomes. Donc c'est quand même plutôt lent. Je ne sais pas si GAMESS utilise le GPU des RPi pour ses calculs ?

Merci Luthaf pour le point de comparaison.

Je vais relancer le calcul demain sur GAMESS avec mon ordi au boulot (un simple PC de bureau). Malheureusement je n'ai pas d'autres points de comparaison. Concernant le GPU, faut que voie, je ne sais pas.

Je vais essayer un exemple "pratique" et je vous tiens au courant des performances.

Dans tout les cas, si j'ai un jour le temps et l'énergie, j'écrirai un tuto de chimie quantique, et vu que GAMESS est gratuit (même si sous licence), je pense proposer de s'en servir pour faire deux ou trois trucs drôle.

Un autre code, open-source pour sa part : http://www.quantum-espresso.org/. Ça fonctionne plutôt pas mal, j'ai un petit projet de métadynamique ab initio avec et ça tourne bien.

Pas mal, je connaissais pas Ceci dit, ça m'as l'air vachement orienté DFT, et si j'écris un truc, je serai obligé de commencer par Hatree-Fock. Donc je note

Si Lapack/BLAS le fait, oui. Mais ça m'étonnerai.

Juste pour donner une comparaison sur un noeud de calcul, pour effectuer la même optimisation :

1

TOTAL WALL CLOCK TIME=       45.7 SECONDS, CPU UTILIZATION IS  99.56%

(Intel sandy bridge, 2.2Ghz¹, 1 coeur², 1Gio de mémoire³)

À noter que c'est la version 2012 de GAMESS, du coup j'ai été obligé d'employer cet input-ci:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

 $CONTRL DFTTYP=M06-L RUNTYP=OPTIMIZE MAXIT=30 MULT=1 $END
 $SYSTEM TIMLIM=525600 MEMORY=65536000 $END
 $BASIS GBASIS=DZV $END
 $SCF DIRSCF=.TRUE. $END
 $STATPT OPTTOL=0.0005 NSTEP=50 $END
 $DATA 
METHANE
C1
C     6.0    -5.29000    -2.22000     0.00000
H     1.0    -4.22000    -2.22000     0.00000
H     1.0    -5.65000    -2.71000     0.88000
H     1.0    -5.65000    -2.71000    -0.86000
H     1.0    -5.65000    -1.21000    -0.01000
 $END

plus aurait été du gachis

C'est pas parallélisé comme code ?

D'après la page officielle, si, bien entendu, pour la plupart des codes (tableau au milieu de la page, chercher les "p"). Quand je dit "plus ça aurait été du gâchis", je veux bien entendu dire que je sais pertinemment que ce genre de calcul est tellement court que ça ne sert à rien que je demande 4 coeurs pour le plaisir, d'autant que je suis soumis à un système de priorité et que moins demande mes jobs, plus vite ils passent

TOTAL WALL CLOCK TIME= 45.7 SECONDS, CPU UTILIZATION IS 99.56% ```

Donc 30x plus rapide que le RPi! Aïe.

Qu'est-ce que tu as dû changer? MULT n'est pas par défaut? Tes valeurs dans $SYSTEM sont calibrées pour ta machine?

Objectivement, j'ai pas trop compris : si je prenais ton input tel quel, il me faisait un calcul d'énergie sans faire d'optimisation. Du coup, j'ai repris m'input, je l'ai repassé dans MacMolPlot, je l'ai ressorti, je l'ai relancé. Les seules différences que je vois, c'est la définition explicite de la multiplicité (et à la limite, d'un peu plus de mémoire, mais il a pas employé tout ça, donc soit). À mon avis, effectivement, sur la version de 2012, MULT est pas par défaut. Par contre, les valeurs de $SYSTEM, comme je te dis, je peux encore mettre ce que je veux, rien de spécifique là derrière

Ok, bizarre…

Sinon sais-tu comment inclure un set de base (par exemple téléchargé ici) directement dans le fichier d'input? Il faut normalement le mettre dans $data mais je n'ai pas trouvé la bonne syntaxe…

Oui (mais je viens de passer 20 minutes à me souvenir comment)

Bon, partons du principe que tu ne parle pas d'effective core potential, avec lequel je t'avoue n'avoir jamais joué (je connais quelqu'un qu'il l'as fait, mais avec Gaussian, donc bon). En fait, c'est tout bête: tu télécharges ta base d'orbitale pour les atomes qui t'intéressent en sélectionnant bien le format "GAMESS (US)", et ensuite tu défini à la suite de chaque atome les primitives que tu met et comment tu les définis:

1
2
3
4
5
6
7

$DATA
ATOME Z.0 0.0 0.0 0.0
[definitions des OA]

[...]

$END

Donc par exemple, pour du H₂S avec la base MIDI (je sais absolument pas à quoi elle correspond, je trouvais juste le nom sympa et j'ai vu qu'elle était pas trop grosse), ça donne :

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44

 $CONTRL SCFTYP=RHF RUNTYP=ENERGY $END

 $DATA 
Energy H2S (MIDI)
C1
SULFUR    16.0    -1.86589     1.59826    -0.00000
S   3
  1   1201.4584000              0.0657650        
  2    181.3921200              0.3839480        
  3     39.4047950              0.6743720        
S   3
  1     52.1390300             -0.0942320        
  2      4.5287990              0.6354680        
  3      1.7549380              0.4315060        
S   2
  1      2.8404370              0.1900420        
  2      0.3814330             -0.6855270        
S   1
  1      0.1387860              1.0000000        
P   3
  1     54.6440710              0.1036730        
  2     12.1229020              0.4581900        
  3      3.2065040              0.6134000        
P   2
  1      0.8632740              0.2294360        
  2      0.1086790              0.3537000        
P   1
  1      0.3181300              1.0000000   

HYDROGEN     1.0    -0.50589     1.59826    -0.00000
S   2
  1      4.5018000              0.0704520        
  2      0.6814440              0.4078260        
S   1
  1      0.1513980              1.0000000   

HYDROGEN     1.0    -2.31922     0.84984     1.04114
S   2
  1      4.5018000              0.0704520        
  2      0.6814440              0.4078260        
S   1
  1      0.1513980              1.0000000 

 $END

Le truc important, c'est de pas donner de ligne $BASIS, et de bien passer une ligne entre chaque définition d'atome, sinon GAMESS te crie dessus comme un malpropre (ce qu'il a tendance à faire souvent, puisqu'il écrit systématiquement en majuscule). Le gros défaut, c'est que tu dois répéter la même chose pour chaque atome. D'un autre coté , il y a ça, que j'ai pas essayé (la flemme d'aller tripoter mes variables d’environnement et le script rungms qui est juste infâme à configurer).

… Du coup, je me suis pris un trip et j'ai tenté le H₂S avec aug-cc-pVQZ. Pour voir

Aaah, je n'avais pas compris qu'il fallait remettre le set pour chaque atome! Merci beaucoup. Du coup on a le même format que le fichier .dat que GAMESS construit en début de run! Merci encore.

C'est la solution que je privilégiais jusqu'à présent. Pour les sets qu'on utilise fréquemment c'est le top! Au final je n'ai changé que la variable EXTBAS dans gms-files.sh (ce n'est pas comme ça qu'on est sensé faire mais ça marche si on donne le nom du fichier complet à EXTBAS)…

Concernant les ECP, prenons par exemple la molécule imaginaire RuCH et le set de base Def2-SVP sur bse.pnl.gov. On constate que seul le Ru possède un ECP. Si on veut l'utiliser, il faut mettre PP=read dans le bloc \$contrl … \$end puis créer un bloc \$ecp … \$end dans le fichier d'input:

Attention: dans \$ecp … \$end, il faut lister tous les atomes dans le même ordre que dans \$data … \$end en spécifiant NONE pour ceux qui n'ont pas d'ECP:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31

 $data
C1
Molecule imaginaire RuCH
C   6.0   ...   ...   ...
H   1.0   ...   ...   ...
Ru  44.0  ...   ...   ...
 $end

 $ecp
  C-ECP NONE
  H-ECP NONE
  RU-ECP GEN     28     3
  2      ----- f-ul potential     -----
    -28.34061627       2     11.3600000        
     -4.94462923       2      5.6800000        
  4      ----- s-ul potential     -----
    209.82297122       2     11.1052693        
     30.65472642       2      5.4147454        
     28.34061627       2     11.3600000        
      4.94462923       2      5.6800000        
  4      ----- p-ul potential     -----
    146.33618228       2      9.7712707        
     24.12787723       2      5.0739908        
     28.34061627       2     11.3600000        
      4.94462923       2      5.6800000        
  4      ----- d-ul potential     -----
     67.51589667       2      7.6714231        
      9.87010415       2      4.1365647        
     28.34061627       2     11.3600000        
      4.94462923       2      5.6800000        
  $end

Bonne chance avec aug-cc-pVQZ, moi ma molécule organique de 54 atomes en B3LYP/6-331+G(d,p) lancée hier soir sur RPi tourne toujours…

Edit: correction de ma réponse concernant l'utilisation d'un ECP

Aaah, je n'avais pas compris qu'il fallait remettre le set pour chaque atome! Merci beaucoup. Du coup on a le même format que le fichier .dat que GAMESS construit en début de run! Merci encore.

C'est ça

Concernant les ECP, ce n'est pas très compliqué: si ton set de base possède un ECP (par exemple Def2-SVP pour le ruthénium tiré de bse.pnl.gov), tu mets cette partie dans le fichier d'input (avant ou après \$data):

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21

$ECP RU-ECP GEN 28 3 2 ----- f-ul potential ----- -28.34061627 2 11.3600000 -4.94462923 2 5.6800000 4 ----- s-ul potential ----- 209.82297122 2 11.1052693 30.65472642 2 5.4147454 28.34061627 2 11.3600000 4.94462923 2 5.6800000 4 ----- p-ul potential ----- 146.33618228 2 9.7712707 24.12787723 2 5.0739908 28.34061627 2 11.3600000 4.94462923 2 5.6800000 4 ----- d-ul potential ----- 67.51589667 2 7.6714231 9.87010415 2 4.1365647 28.34061627 2 11.3600000 4.94462923 2 5.6800000 $END

et tu ajoutes PP=READ dans \$contrl … $end.

J'ai pas de raison d'essayer dans l'absolu, mais merci

Oh la grosse b…ase d'orbitale atomique. C'est pas commun, ça, tu cherches un truc en particulier ?

J'essaie juste de reproduire une optimisation de géométrie trouvée dans un article d'organique. C'est pour voir si RPi/GAMESS est suffisant pour produire quelque chose d'"utile" dans un temps raisonnable.

Ah, oki Bah tient nous au courant.

(moi, avec l'habitude, je dis 4 jours).

J'en suis à 7 itérations en 48h… Du coup ton estimation sera proche!

Bah j'ai bossé pendant 1 an avec des systèmes qui tournaient autour de 50 atomes, donc j'ai l'habitude de jouer avec ce genre de chose (ceci dit, c'était pas avec GAMESS, donc faut voir l'algo qui est derrière qui est peut être pas optimisé pareil).

Je viens d'écrire un petit script python qui transforme

1
2
3
4
5

C    6.0    -5.29    -2.22    0.00
H    1.0    -4.22    -2.22    0.00
H    1.0    -5.65    -2.71    0.88
H    1.0    -5.65    -2.71    -0.86
H    1.0    -5.65    -1.21    -0.01

en

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32

$data
...
C1
C    6.0    -5.29    -2.22    0.00
S   2
1     27.38503303             0.430128498
2      4.874522052            0.678913531
L   2
1      1.136748189            0.049471769            0.511540708
2      0.288309360            0.963782408            0.612819896

H    1.0    -4.22    -2.22    0.00
S   2
1      1.309756377            0.430128498
2      0.233135974            0.678913531

H    1.0    -5.65    -2.71    0.88
S   2
1      1.309756377            0.430128498
2      0.233135974            0.678913531

H    1.0    -5.65    -2.71    -0.86
S   2
1      1.309756377            0.430128498
2      0.233135974            0.678913531

H    1.0    -5.65    -1.21    -0.01
S   2
1      1.309756377            0.430128498
2      0.233135974            0.678913531

$end

à partir de

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62

!  STO-2G  EMSL  Basis Set Exchange Library   13/12/14 12:26
! Elements                             References
! --------                             ----------
!  H - He: W.J. Hehre, R.F. Stewart and J.A. Pople, J. Chem. Phys. 2657 (1969).
! Li - Ne: 
! Na - Ar: W.J. Hehre, R. Ditchfield, R.F. Stewart, J.A. Pople, J. Chem. Phys.
!  K - Kr: 52, 2769 (1970).
!   


$DATA
HYDROGEN
S   2
 1      1.309756377            0.430128498      
 2      0.233135974            0.678913531      

CARBON
S   2
 1     27.38503303             0.430128498      
 2      4.874522052            0.678913531      
L   2
 1      1.136748189            0.049471769            0.511540708      
 2      0.288309360            0.963782408            0.612819896      

NITROGEN
S   2
 1     37.8964750              0.4301280        
 2      6.7455540              0.6789140        
L   2
 1      1.4610890              0.0494720              0.5115410        
 2      0.3705710              0.9637820              0.6128200        

OXYGEN
S   2
 1     49.9809710              0.4301280        
 2      8.8965880              0.6789140        
L   2
 1      1.9452370              0.0494720              0.5115410        
 2      0.4933630              0.9637820              0.6128200        

PHOSPHOROUS
S   2
 1    179.0948740              0.4301280        
 2     31.8787970              0.6789140        
L   2
 1     10.8341890              0.0494720              0.5115410        
 2      2.7478370              0.9637820              0.6128200        
L   2
 1      0.5872910             -0.2983990              0.3480470        
 2      0.2014940              1.2279830              0.7222520        

CHLORINE
S   2
 1    229.9441040              0.4301280        
 2     40.9299350              0.6789140        
L   2
 1     15.0576100              0.0494720              0.5115410        
 2      3.8190080              0.9637820              0.6128200        
L   2
 1      0.8883460             -0.2983990              0.3480470        
 2      0.3047830              1.2279830              0.7222520        
$END

Je vais l'améliorer un peu et le partager.