Faire Vectoriser une boucle avec icpc

**reptils** · 08/03/2010, 15h15

Bonjour,
Voulant mimiquer le comportement de Fortran, je suis en train de faire mumuse en C++ pour faire une classe Array (d'abord en 1d).
Je pensais "naivement" que le compilateur intel pourrait facilement m'aider à vectoriser une simple boucle où j'utilise un opérateur () pour accéder à un champ:

bout du code de la class Array_1d:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

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
template<class T,int Offset=DEF_OFFSET>
class Array_1d
{
    private:
        uint64_t n1;
        int64_t o1;
 
        T* RESTRICT ptr;
 
        // Forbid the "pass by value" in function calls
        Array_1d(const Array_1d<T,Offset>& rhs);
    public:
        // Constructors
        Array_1d(const uint64_t d1=0);
        Array_1d(const Range& r1);
 
        // Destructor
        ~Array_1d();
 
        // Array_1d<T> operators
        const Array_1d<T,Offset>& operator=(const Array_1d<T,Offset>& rhs);
        const Array_1d<T,Offset>& operator=(const T value);
 
        inline const T& operator()(const int64_t& i) const
        {
#ifdef CHECKBOUND
            if (i<o1 || i>int64_t(o1+n1-1))
                throw Error(1,"Out of bound");
#endif
            return ptr[i-o1];
        }
        // OPERATEUR QUI DOIT ETRE VECTORISER EN BOUCLE SELON MOI !!!
        inline T& operator()(const int64_t& i)
        {
#ifdef CHECKBOUND
            if (i<o1 || i>int64_t(o1+n1-1))
                throw Error(1,"Out of bound");
#endif
            return ptr[i-o1];
        }
 
        // Methods
        inline const T* RESTRICT get_pointer() const {return ptr;}
        inline T* RESTRICT get_pointer() {return ptr;}
 
        inline uint64_t size1() const{return n1;}
 
        inline Range range1() const{return Range(o1,o1+n1-1);}
 
        void resize(const uint64_t n_d1);
        void resize(const Range& r1);
        void free();
};

RESTRICT est définit en tant que restrict avec le compilo intel, et je compile avec -O3 -xSSSE3 -restrict.

En ce qui concerne une simple boucle d'appel à l'opérateur, je fais:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
 
    int64_t i;
    int64_t min,max;
    Array_1d<float> array;
 
    array.resize(Range(-200,400));
    min=array.range1().min();
    max=array.range1().max();
 
// Ne vectorise pas ???    
    for(i=min;i<=max;++i)
        array(i)=1.5;
 
 
    float* ptr=array.get_pointer();
 
// Vectorise !!!
    for(i=min;i<=max;++i)
        ptr[i]=1.5;

Est ce que vous avez des idées? ou est ce qu'on arrive dans les limites du/des compilos C++ qui s'arretent juste à l'inline sans passer à l'étape d'optimisation d'après (la vectorisation)?

Faisant pas mal de calcul scientifique je trouve dommage de ne pas pouvoir vectoriser ce genre de chose... pour les performances c'est très important.

à vos commentaires !!!

**Matthieu Brucher** · 09/03/2010, 08h52

Salut,

Peux-tu poster un code complet compilable ? Il manque des éléments pour l'instant.
Si tu veux faire mumuse, je te conseille de regarder Eigen qui est la référence actuelle en terme de calcul matriciel en C++ (et même limite en terme de calcul matriciel tout court car les perfs sont équivalentes à la MKL et GOTO BLAS, même en multithreadé).

**reptils** · 09/03/2010, 13h44

Je ne connaissais pas Eigen... mais je me méfie des annonces de perf battant du BLAS etc. On peut toujours en interne gérer des calculs et les rendre aussi performant que le C ou le Fortran.
Mais c'est moins le cas dans ce que je montre dans mon exemple. Là c'est offrir un accès à une matrice comme si on faisait du Fortran ou accéder à un tableau C directement. J'ai regardé le code source de Blitz++ ou de tvmet et ils font quasiment la même chose que ce que j'ai présenté... mais ça ne vectorise pas en externe.

Je fournirai un code source qui compile dans la journée. Mais je reste sceptique sur l'idée d'arriver à faire vectoriser ce machin

**Goten** · 09/03/2010, 13h46

Blitz++ (qui n'est plus maintenu) utilise des expressions templates principalement ... Me rappelle plus si y'a de la vectorisation.

Eigen par contre fait les deux.

**Matthieu Brucher** · 09/03/2010, 15h53

Envoyé par reptils

Je ne connaissais pas Eigen... mais je me méfie des annonces de perf battant du BLAS etc. On peut toujours en interne gérer des calculs et les rendre aussi performant que le C ou le Fortran.

http://eigen.tuxfamily.org/index.php?title=Benchmark
J'ai vu sur la ML de dev d'Eigen des benchs plus à jour.

Envoyé par reptils

Mais c'est moins le cas dans ce que je montre dans mon exemple. Là c'est offrir un accès à une matrice comme si on faisait du Fortran ou accéder à un tableau C directement. J'ai regardé le code source de Blitz++ ou de tvmet et ils font quasiment la même chose que ce que j'ai présenté... mais ça ne vectorise pas en externe.

C'est faisable sans aucun problème avec vectorisation à la clé (sachant que dans le cas d'une addition, la vectorisation n'a aucun impact

)

Envoyé par reptils

Je fournirai un code source qui compile dans la journée. Mais je reste sceptique sur l'idée d'arriver à faire vectoriser ce machin

On va y arriver. Il suffit d'analyser ce qu'ICPC sort avec -vec-report=2, et on aura déjà une bonne idée de pourquoi il ne vectorise pour l'instant pas. Mais il faut un code simple qui compile

(pff, tu es en vacances, tu devrais te reposer plutôt que de faire du C++, à moins de vouloir convertir HC ?

)

**reptils** · 09/03/2010, 16h58

Envoyé par Matthieu Brucher

http://eigen.tuxfamily.org/index.php?title=Benchmark
J'ai vu sur la ML de dev d'Eigen des benchs plus à jour.

Oui j'ai vu mais comme je te le dis, ce sont des paternes qui sont gérées en interne

C'est faisable sans aucun problème avec vectorisation à la clé (sachant que dans le cas d'une addition, la vectorisation n'a aucun impact

)

???

On va y arriver. Il suffit d'analyser ce qu'ICPC sort avec -vec-report=2, et on aura déjà une bonne idée de pourquoi il ne vectorise pour l'instant pas. Mais il faut un code simple qui compile

La log peut être analysée avec le tar que je donne ci dessous... il me sort depuis hier soir la même chose (dependency FLOW ANTI... blablabla), pas moyen de faire autrement. Pour le voir il faut monter le niveau du vec-report à 3.

(pff, tu es en vacances, tu devrais te reposer plutôt que de faire du C++, à moins de vouloir convertir HC ?

)

Non c'est juste pour faire mumuse, je voulais faire des tests par la suite...

Tu as l'archive de mon test ici:
http://reptils.free.fr/Developpez.tgz

La compile me donne:

Processing main.cpp
icpc -c -O3 -xssse3 -g -vec-report2 main.cpp -DDEF_OFFSET=0 -DRESTRICT=restrict -restrict
main.cpp(18): (col. 5) remark: loop was not vectorized: existence of vector dependence.
main.cpp(24): (col. 5) remark: LOOP WAS VECTORIZED.
main.cpp(28): (col. 10) remark: LOOP WAS VECTORIZED.
icpc -O3 -xssse3 -g -vec-report2 -o test.bin main.o

la ligne 18 est celle qui à mon avis devrait être vectorisée si l'inline est bien fait???

EDIT: j'oubliais, j'ai regarder le code source de Eigen, je n'ai pas l'impression que ça fasse ce que je décrit.

**Matthieu Brucher** · 09/03/2010, 20h17

Envoyé par reptils

Oui j'ai vu mais comme je te le dis, ce sont des paternes qui sont gérées en interne

Pas compris ce que tu voulais dire...

Envoyé par reptils

???

Pour une addition, ton facteur limitant, c'est la bande passante mémoire. Pour une multiplication, c'est différent.

Envoyé par reptils

La log peut être analysée avec le tar que je donne ci dessous... il me sort depuis hier soir la même chose (dependency FLOW ANTI... blablabla), pas moyen de faire autrement. Pour le voir il faut monter le niveau du vec-report à 3.

2 ne suffit pas ?
En tout cas, ce sont ces éléments qui indiquent ce qu'il faut faire.

Envoyé par reptils

la ligne 18 est celle qui à mon avis devrait être vectorisée si l'inline est bien fait???

Pas forcément. Il y a des histoires d'assumed dependencies et proved dependencies. Si c'est de l'assumed, un #pragma ivdep résoud le pb, normalement.

Envoyé par reptils

EDIT: j'oubliais, j'ai regarder le code source de Eigen, je n'ai pas l'impression que ça fasse ce que je décrit.

En meêm temps, je n'ai pas tout à fait saisi ce que tu voulais faire

Quand tu dis comme en Fortran, c'est en pouvant travailler sur [-n,n] par exemple au lieu de [0,2n] ?

**Matthieu Brucher** · 09/03/2010, 21h02

Je commence à regarder, mais sous Windows. Tu as 2-3 petits problèmes, comme l'utilisation de stdint.h qui n'est pas standard (il faut passer par Boost normalement), ou l'utilisation des vieux en-têtes. Il faut que je vois encore comment je génère les rapports de vectorisation sous Windows...

**reptils** · 09/03/2010, 21h13

Ce que je dis c'est que comme pour la ligne qui vectorise dans mon cas

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

array=1.5

C'est une gestion interne à ma classe quiu fait que cette instruction vectorise. Et c'est le cas pour l'ensemble des instructions de Eigen. C'est en ça que je dis que c'est facile de faire vectoriser des choses.

Ensuite, je sais pour l'addition, mais bon, toute vectorisation est bonne à prendre

tu ne sais pas par exemple si tes données sont déjà en cache, dans ce cas tu es content d'avoir tes instructions qui sont vectorisées.

Pas forcément. Il y a des histoires d'assumed dependencies et proved dependencies. Si c'est de l'assumed, un #pragma ivdep résoud le pb, normalement.

Déjà testé, mais n'a pas marché.

En meêm temps, je n'ai pas tout à fait saisi ce que tu voulais faire
Quand tu dis comme en Fortran, c'est en pouvant travailler sur [-n,n] par exemple au lieu de [0,2n] ?

Oui, je veux pouvoir définir la "shape" que je veux. Faire des trucs du genre:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
 
for(int i=xmin,i<=xmax;++i)
{
// do something with array(i)
}

avec xmin et xmax qui valent n'importe quoi, positif, négatif... tout en pouvant faire mes instructions se vectoriser

Non mais je pense que c'est une limitation des compilos qui n'arrivent pas à avoir une visibilité nécessaire pour vectoriser, malgré l'inline.
C'est dommage que le c++ n'offre pas ce genre de souplesse quant meme...

**Matthieu Brucher** · 09/03/2010, 21h49

Envoyé par reptils

Ce que je dis c'est que comme pour la ligne qui vectorise dans mon cas

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

array=1.5

C'est une gestion interne à ma classe quiu fait que cette instruction vectorise. Et c'est le cas pour l'ensemble des instructions de Eigen. C'est en ça que je dis que c'est facile de faire vectoriser des choses.

Je dois t'avouer que sur ma version du compilateur, il ne comprend pas la version interne à la classe et n'arrive pas à la vectoriser... En revanche, il comprend maintenant très bien la boucle ligne 23-24 et la vectorise.
Pour la première et la troisième vectorisation, l'assumed dependency levé, il reste une boucle trop complexe pour la version que j'ai

Je pense que sur une version plus récente, ça passerait, il faut que je teste.

Envoyé par reptils

Ensuite, je sais pour l'addition, mais bon, toute vectorisation est bonne à prendre

tu ne sais pas par exemple si tes données sont déjà en cache, dans ce cas tu es content d'avoir tes instructions qui sont vectorisées.

On est d'accord

Envoyé par reptils

Déjà testé, mais n'a pas marché.

Bizarre, sur une version plus ancienne, ça fonctionne et il y a qqs mois, un patch a été ajouté au compilo pour améliorer encore la qualité des messages. Si tu as mis un tel pragma, au moins tu devrais voir une meilleure raison que celle affichée.

Envoyé par reptils

Oui, je veux pouvoir définir la "shape" que je veux. Faire des trucs du genre:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
 
for(int i=xmin,i<=xmax;++i)
{
// do something with array(i)
}

avec xmin et xmax qui valent n'importe quoi, positif, négatif... tout en pouvant faire mes instructions se vectoriser

Non mais je pense que c'est une limitation des compilos qui n'arrivent pas à avoir une visibilité nécessaire pour vectoriser, malgré l'inline.

En même temps, ton code a des petits problèmes tout de même.
- le mot-clé inline est inutile dans une classe template si tu définis le code à l'intérieur de la classe, ce que tu fais avec (), mais pas pour les autres fonctions qui risquent de se retrouver à plusieurs instanciations

- Ta classe est template avec l'offset comme paramètre template, tu ne devrais donc pas pouvoir changer l'offset de la classe, or tu le peux dans ton cas. il est possible que le compilateur ait du mal à simplifier cette opération. En fait, ton paramètre Offset est en réalité complètement inutile

Envoyé par reptils

C'est dommage que le c++ n'offre pas ce genre de souplesse quant meme...

Il est certain que les compilateurs C et C++ ne permettent pas encore le même niveau que Fortran, mais on est tout de même très proche sur les bibliothèques de référence, aux ranges() près (range() qui ne sont pas non plus très souvent utilisés en Fortran...).

Update : il semblerait que ma version est la dernière disponible sous Windows (j'ai même une version plus récente que celle proposée par Intel sur son site

)

**Joel F** · 09/03/2010, 21h54

La plupart des algorithmes de vectorisations des compilateurs ont du mal avec des bornes de boucles non constantes ou non inférables comme constantes ainsi que pour des nids de boucles multi-dimensionels.

**Matthieu Brucher** · 09/03/2010, 22h01

Exact.

D'ailleurs, au lieu de o1, je préfèrerai avoir un ptr_offset qui vaut ptr - Offset. Dans ce cas, si la troisième boucle est vectorisée, la première devrait l'être si on passe par ce pseudo pointeur (malheureusement, je ne peux pas vérifier sur le coup

(()

**Joel F** · 09/03/2010, 22h04

Après, c'est plus productif de sortir son couteau suisse SSEx ou AltiVec est écrire les 2-3 lignes nécessaires plutôt que de rendre le code illisible pour essayer que le compilo s'y retrouve.

**Matthieu Brucher** · 09/03/2010, 22h10

Pas tout à fait d'accord

Ecrire en SSEx, c'est plus difficile à lire pour le commun des mortels que d'ajouter un restrict et un pragma...

Avec la version sous Windows, il me faut un pointeur restrict temporaire pour que le =(value) soit vectorisé. Comme il est impossible de faire cela avec le premier, je suis coincé à cause de ma version. Mais je testerai demain sous Linux, ça devrait mieux passer. Tu as essayé avec quelle version du compilo ?

**Joel F** · 09/03/2010, 22h17

Envoyé par Matthieu Brucher

Pas tout à fait d'accord

Ecrire en SSEx, c'est plus difficile à lire pour le commun des mortels que d'ajouter un restrict et un pragma...

pour une certaine définition de difficile :o

**reptils** · 10/03/2010, 20h46

Pour info ma version de compilo est 11.1.064/Linux x86_64 donc pas très vieille

Envoyé par Matthieu Brucher

Ta classe est template avec l'offset comme paramètre template, tu ne devrais donc pas pouvoir changer l'offset de la classe, or tu le peux dans ton cas. il est possible que le compilateur ait du mal à simplifier cette opération. En fait, ton paramètre Offset est en réalité complètement inutile

Non, c'est dans le cas où je souhaite faire un resize(N) là Offset me sert à me retrouver dans une configuration C (démarrage en 0) ou Fortran (démarrage en 1).

Comme je le disais au début du topic j'essaie de mimiquer au mieux les comportement du Fortran (avec en prime la possibilité d'indexer à partir de 0 c'est user define).

Envoyé par Matthieu Brucher

Il est certain que les compilateurs C et C++ ne permettent pas encore le même niveau que Fortran, mais on est tout de même très proche sur les bibliothèques de référence, aux ranges() près (range() qui ne sont pas non plus très souvent utilisés en Fortran...).

Faux! C'est très pratique et utilisé.
Exemple simple, une décomposition en sous-domaine si tu fais du MPI, chaque MPI peut avoir son tableau de sous-domaine avec une indexation globale (un domaine d'application que tu connais très bien...)

Mais de toute façon ce n'est pas l'offset qui pose le problème de la vectorisation (encore heureux!) Le compilo gère très bien le décalage d'offset.
Visiblement ce qu'il n'aime pas, c'est de ne pas avoir assez d'infos pour s'assurer que la classe ne bouge pas entre temps au cours du déroulage de boucle.

Quoi qu'il en soit je pense que des efforts doivent être fait à ce niveau là (je parle des compilos). Parce que les librairies de référence, dont tu fais référence

offrent un subset de routines pour un scope donnée. Dans le cas simple que j'essaie de mettre en place c'est la vectorisation d'une boucle qui est simplexe, et qui peut être utilisée dans n'importe quel contexte. Certaines applications numériques n'utilisent pas de MKL par exemple, et là tu es bien content de pouvoir quant même faire du code qui vectorise.

Envoyé par Joel F

La plupart des algorithmes de vectorisations des compilateurs ont du mal avec des bornes de boucles non constantes ou non inférables comme constantes ainsi que pour des nids de boucles multi-dimensionels.

Ben dans mon cas j'ai sorti les bornes pour les mettre dans des variables min et max. Donc le problème ne vient pas de là. Surtout que si je fais directement (en mettant les membres public):

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
for(i=min;i<=max;++i)
    array.ptr[i-array.o1]=1.5;

ça se vectorise.

Au final, on peut faire en sorte que le C++ (écrit plutôt en mode C) vectorise mais il faut être très prudent à comment l'utiliser. Je voulais juste, dans mon cas, arriver à créer un ensemble de classes qui donnent la même souplesse de syntaxe que le fortran pour l'accès aux tableau 1d,2d,3d... dans le calcul scientifique. Mais bon, peut être qu'il faut attendre que les compilos mûrissent ou que la norme du C++ apporte quelque chose à ce niveau. Malheureusement le calcul scientifique n'est pas trop une priorité du comité c++

**Matthieu Brucher** · 10/03/2010, 21h01

Envoyé par reptils

Pour info ma version de compilo est 11.1.064/Linux x86_64 donc pas très vieille

Il est possible que la 69 soit plus expressive tout de même.

Envoyé par reptils

Non, c'est dans le cas où je souhaite faire un resize(N) là Offset me sert à me retrouver dans une configuration C (démarrage en 0) ou Fortran (démarrage en 1).

Comme je le disais au début du topic j'essaie de mimiquer au mieux les comportement du Fortran (avec en prime la possibilité d'indexer à partir de 0 c'est user define).?

Justement, tu ne devrais pas avoir ce paramètre en template. C'est une valeur par défaut qui est overridée dès le départ. Avec ta macro, tu fais tout aussi bien, mais en plus clair pour l'utilisateur.
Les classes matricielles qui utilisent des variables template sont assez claires : leur contenu a un impact sur toute la durée de vie de l'objet, pas juste à la construction. Une valeur par défaut dans le constructeur est largement suffisante pour ton objectif.

Envoyé par reptils

Faux! C'est très pratique et utilisé.
Exemple simple, une décomposition en sous-domaine si tu fais du MPI, chaque MPI peut avoir son tableau de sous-domaine avec une indexation globale (un domaine d'application que tu connais très bien...)

Pas vraiment... Je te rappelle que justement ce genre de code est plus souvent en simili Fortran 77 qu'en Fortran 90 avec ce que tu dis. De mon côté, je préfèrerai aussi que ce soit plus utilisé, je trouve ça clair et élégant (et comme tu le sous-entends, on sait bien mieux ce qu'on fait

)

Envoyé par reptils

Mais de toute façon ce n'est pas l'offset qui pose le problème de la vectorisation (encore heureux!) Le compilo gère très bien le décalage d'offset.
Visiblement ce qu'il n'aime pas, c'est de ne pas avoir assez d'infos pour s'assurer que la classe ne bouge pas entre temps au cours du déroulage de boucle.

Les règles du C++ fait que le compilateur sait que c'est possible. Malheureusement, je n'ai pas eu le temps aujourd'hui de tester la version que j'ai mise en ligne. De manière globale, il faut simplifier la vie du compilateur. Tu as essayé en remplaçant o1 par un ptr_offset = ptr - offset ?

Envoyé par reptils

Quoi qu'il en soit je pense que des efforts doivent être fait à ce niveau là (je parle des compilos). Parce que les librairies de référence, dont tu fais référence

offrent un subset de routines pour un scope donnée. Dans le cas simple que j'essaie de mettre en place c'est la vectorisation d'une boucle qui est simplexe, et qui peut être utilisée dans n'importe quel contexte. Certaines applications numériques n'utilisent pas de MKL par exemple, et là tu es bien content de pouvoir quant même faire du code qui vectorise.

Ce qui est vraiment bizarre, c'est que ça passe pour Eigen sans problème, il faut vraiment que je regarde ce qui se passe

Envoyé par reptils

Ben dans mon cas j'ai sorti les bornes pour les mettre dans des variables min et max. Donc le problème ne vient pas de là. Surtout que si je fais directement (en mettant les membres public):

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
for(i=min;i<=max;++i)
    array.ptr[i-array.o1]=1.5;

ça se vectorise.

Au final, on peut faire en sorte que le C++ (écrit plutôt en mode C) vectorise mais il faut être très prudent à comment l'utiliser. Je voulais juste, dans mon cas, arriver à créer un ensemble de classes qui donnent la même souplesse de syntaxe que le fortran pour l'accès aux tableau 1d,2d,3d... dans le calcul scientifique. Mais bon, peut être qu'il faut attendre que les compilos mûrissent ou que la norme du C++ apporte quelque chose à ce niveau. Malheureusement le calcul scientifique n'est pas trop une priorité du comité c++

C'est aux compilateurs de mieux se débrouiller. Je ne pense pas que le C ou le C++ offre un jour une écriture semblable à Fortran/Matlab/Python.
C'est vraiment bizarre, parce qu'ils arrivent à vectoriser des codes plus complexes

Comme pour tout, il faut trouver le petit déclic pour savoir comment trouver la solution (comme lorsqu'on vectorise un stencil en Fortran, on le découpe en sous-stencil vectorisable) et après, on peut l'appliquer partout

**Matthieu Brucher** · 11/03/2010, 09h03

Bon, je n'ai pas réussi à lui faire comprendre ce que je voulais pour l'instant...

De toute manière, il y a un point important à noter. Si on veut de la perf en Fortran, on passe par l'écriture par array. A ce moment, il n'y a qu'une seule opération faite sur la totalité des données (par ex Z = aX + bY). Pour faire de même en C++ sans variables intermédiaires, il faut aussi passer par des expression templates qui sont faites automatiquement avec Eigen et autres Blitz++.

Si je suis d'accord sur le fait que les compilateurs devraient pouvoir vectoriser les codes simples comme le premier, ils sont au moins capables de le faire avec l'écriture tableau = valeur. Pour des opérations plus complexes non disponibles avec une bibliothèque dédiée, il n'est pas certain que Fortran soit capable de faire mieux.

Tout cela pour dire qu'en modifiant Eigen (ajouter un argument offset), il est possible de réaliser ce que tu souhaites (AMHA).

**JolyLoic** · 11/03/2010, 10h26

Envoyé par reptils

Malheureusement le calcul scientifique n'est pas trop une priorité du comité c++

Détrompes toi... Après tous, c'est plus dans les domaines scientifiques que par exemple dans la gestion qu'on a besoin des perfs que peut apporter le C++.

Pas mal de contributeurs au standard sont liés à des universités ou à des labos (fermilab par exemple). Ils n'ont peut-être pas le même poids que les fabricants de compilateurs ou de bibliothèques standard, mais parmi les représentants des utilisateurs, ce sont probablement les plus nombreux.

Et depuis peu, il y a au sein du comité de normalisation une groupe d'intérêt autour du calcul scientifique qui s'est formé, avec comme but de promouvoir des modifications du langage ou des bibliothèques standard dans ce sens. Par contre, ce groupe étant récent, je ne crois pas qu'il ait eu des impacts concrets.

Mais dans les évolutions du langage, il y a quand même des points autour de ça (fonction mathématiques spéciales, générateur de nombre aléatoires, on a failli avoir calcul par intervalle...).

**reptils** · 11/03/2010, 19h50

Envoyé par Matthieu Brucher

De toute manière, il y a un point important à noter. Si on veut de la perf en Fortran, on passe par l'écriture par array. A ce moment, il n'y a qu'une seule opération faite sur la totalité des données (par ex Z = aX + bY)....

Faux, c'est une grosse idée reçue, et je pense qu'il ne faut surtout pas utiliser la syntaxe array. C'est un truc de feignasses qui pensent que le compilo est super intelligent et que ça permet d'éviter l'écriture de boucles.

En fait on croit que c'est plus simple pour le compilo alors que dans les faits, il semble qu'il a plus de mal à se dépatouiller surtout si tu utilises des tableaux 2d 3d avec des écritures alambiquées... (1d il gère bien, encore heureux).

Si tu veux de la perf il faut expliciter au mieux tes boucles. Plus tu donnes de l'info au compilo et mieux il vectorise.
Sur des algo type tracé de rayon (si tu vois de quoi je parle) ils ont des facteurs 2 à 4 en virant la syntaxe array sur ces portions de code. De mon coté on garde que des expressions simples. (ex: array(:,:,: )=cste)

Envoyé par Matthieu Brucher

Pour des opérations plus complexes non disponibles avec une bibliothèque dédiée, il n'est pas certain que Fortran soit capable de faire mieux.

Je pense plutôt le contraire. La MKL elle même est écrite en C et asm. Mais comme je te le disait au début du topic, c'est une gestion interne, donc les Blitz++, Eigen peuvent faire tout aussi bien (enfin presque).

Le problème c'est que si tes algos sortent du scope de ces librairies, tu aimerais bien pouvoir écrire quant même du code performant, avec des expressions simples:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
 
for(k=zmin;k<=zmax;++k)
{
    for(j=ymin;j<=ymax;++j)
    {
        for(i=xmin;i<=xmax;++i)
        {
            array1(i,j,k)=alpha*array2(i-4,j+3,k+n) OP beta*array3(i+c,j-d,k+z) OP array4(i,j,k)...
         }
    }
}

Tu ne trouveras pas de librarie pour t'offrir forcement ce qu'il y a au dessus. Si le compilo comprend que tu n'es pas aliasé et que ton expression est vectorisable, alors il le fait et tu es gagnant.

Matthieu remarque bien que dans tous les cas de figure je ne parle pas des problèmes de gestion de la mémoire. De mon point de vue c'est pas quelque chose qui est lié au langage. Que ce soit du C/C++/Fortran c'est une question d'organisation de boucle et de la gestion des tableaux. En revanche si ton code est CPU bound, ben tes accélérations tu vas les chercher dans le SIMD.

Envoyé par JolyLoic

Mais dans les évolutions du langage, il y a quand même des points autour de ça (fonction mathématiques spéciales, générateur de nombre aléatoires, on a failli avoir calcul par intervalle...).

Il y a en effet ce point, c'est une chose et c'est très important. Mais je pense qu'on devrait être capable d'écrire du code (comme celui que j'ai présenté) et avoir le compilo capable de résoudre le code pour le vectoriser. Car si un compilo est capable de vectoriser sur les operator() de ce genre de classe (array_1d, array_2d,array_3d,...) alors il n'y aurait quasiment plus lieu de mettre en avant le fortran comme étant plus avantagieux que le C++ pour les perfs et la facilité d'écriture de code pour un scientifique non développeur chevronné.
On aurait alors un langage ultra-complet:
- Possibilité d'utiliser les expression Template (haut niveau)
- Possibilité d'écrire à la Fortran (moyen niveau)
- Possibilité de partir sur un mode plus proche du C (très bas niveau type kernel de calcul)

Faire Vectoriser une boucle avec icpc

C++

Discussions similaires

Partager

Partager