Multiplication efficace des tenseurs

Question

J'ai une matrice qui est une représentation d'un tenseur de dimension supérieure qui pourrait en principe être à N dimensions, mais chaque dimension a la même taille. Disons que je veux calculer ce qui suit :

et C est stocké sous forme de matrice via

où il existe une correspondance entre ij et I et kl et J.

Je peux le faire avec des boucles for imbriquées où chaque dimension de mon tenseur a une taille de 3 via

for (int i=0; i<3; i++){
    for (int j=0; j<3; j++){
        I = map_ij_to_I(i,j);
        for (int k=0; k<3; k++){
            for (int l=0; l<3; l++){
                J = map_kl_to_J(k,l);
                D(I,J) = 0.;
                for (int m=0; m<3; m++){
                    for (int n=0; n<3; n++){
                        M = map_mn_to_M(m,n);
                        D(I,J) += a(i,m)*C(M,J)*b(j,n);
                    }
                }
            }
        }
    }
}

mais c'est plutôt désordonné et pas très efficace. J'utilise la bibliothèque de la matrice d'Eigen et je pense qu'il existe probablement un meilleur moyen de faire cela qu'une boucle for ou le codage de chaque entrée séparément. J'ai essayé la bibliothèque tensorielle non supportée et j'ai trouvé qu'elle était plus lente que mes boucles explicites. Avez-vous des idées ?

En guise de question bonus, comment pourrais-je calculer efficacement quelque chose comme ce qui suit ?

Answer 1

L'optimiseur de votre compilateur fait beaucoup de travail pour vous sous le capot. Pour une fois, les boucles avec un nombre constant d'itérations sont déroulées. C'est peut-être la raison pour laquelle votre code est plus rapide que celui de la bibliothèque.

Je vous suggère de jeter un coup d'œil à l'assemblage produit avec les optimisations activées pour vraiment comprendre où vous pouvez optimiser et à quoi ressemble réellement votre programme une fois compilé.

Ensuite, bien sûr, vous pouvez penser à des implémentations parallèles, soit sur le CPU (plusieurs threads), soit sur le GPU (cuda, OpenCL, OpenAcc, etc.).

Pour ce qui est de la question bonus, si vous envisagez de l'écrire sous la forme de deux boucles imbriquées, je vous suggère de réorganiser l'expression de façon à ce que le terme a_km se trouve entre les deux sommes. Il n'est pas nécessaire d'effectuer cette multiplication à l'intérieur de la somme interne car elle ne dépend pas de n. Bien que cela ne donne probablement qu'un léger avantage en termes de performances sur les CPU modernes...