1
votes

Évaluation du SUM (..) de Fortran sans tableau temporaire?

Au niveau de la lisibilité, je trouve préférable d'écrire

momentum = 0.0d0
do i = 1, size(elements,1)
    momentum = momentum + calculateMomentum(elements(i))
end do

plus de

momentum = sum( [( calculateMomentum(elements(i)), i=1, size(elements,1) )] )

car la première version a la forme de définir la valeur de l' momentum , tandis que la seconde correspond à une instruction d'accumulation de plus bas niveau. La différence devient plus prononcée dans le code du monde réel avec plus d'arguments et éventuellement plusieurs indices.

Cependant, la première version alloue un tableau temporaire. En tant que programmeur humain, je sais que cela pourrait être optimisé, alors je me demandais si Fortran propose une syntaxe, qui permet de calculer la somme sans boucle explicite ni tableau temporaire.


Mise à jour

Pas de telle construction. Il semble qu'il n'y ait pas de syntaxe, comme ce que je demandais. Voir la réponse de Vladimir F.

C'est important, mais moins que je ne le pensais. J'ai fait mon propre benchmark ( pastebin , embarqué ) en utilisant la multiplication matricielle en utilisant plusieurs variantes.

  1. C(i,j) = C(i,j) + A(i,k) * B(k,j) était le plus lent, probablement en raison de l'accès inutile au tableau à chaque étape.

  2. C(i,j) = sum ([(A (i, k) * B (k, j), k = 1, N)]) `était environ 10-20% plus rapide que (1), malgré le tableau temporaire .

  3. tmp = tmp + A(i,k) * B(k,j) , c'est-à-dire en utilisant une variable d'accumulation temporaire, était environ 20% plus rapide que (2).

  4. C = matmul(A,B) était de loin le plus rapide, comparé à (3) par un facteur de 25 pour les matrices 500x500, passant à 50 pour 2000x2000, tandis que la vitesse relative des autres variantes est restée à peu près la même.

Conclusion: lorsque la tâche ne peut pas être exprimée facilement dans des fonctions de bibliothèque ou intrinsèques optimisées, la variante de sum a des performances viables et ne doit être optimisée que si les performances comptent vraiment à un tel degré dans cette partie du code.


0 commentaires

3 Réponses :


2
votes

Le fait qu'un tableau temporaire soit alloué ou net dépend des optimisations du compilateur. L'allocation de pile est de toute façon presque gratuite. La copie des valeurs prendra probablement plus de temps.

Le compilateur peut optimiser les étapes inutiles s'il peut s'assurer que le résultat sera le même. Cependant, il n'y a pas de syntaxe spéciale pour cela. Fortran essaie généralement de rester loin de l'implémentation réelle et laisse beaucoup sur le compilateur.


5 commentaires

Il y a un constructeur de tableau dans l'instruction, mais cela peut être évité.


En Python, il y aurait sum([f(x) for x in xs]) vs sum(f(x) for x in xs) . La première version crée une liste inutile, comme Fortran crée un tableau inutile. J'espérais en quelque sorte m'assurer que l'intrinsèque puisse faire quelque chose de plus similaire à r = 0; do i=1,size(xs); r = r + f(xs(i)); end do tout en conservant le style de code plus déclaratif.


@kdb Il n'y a rien de tel à Fortran. Pas d'itérateurs, pas de fermetures, pas d'évaluation paresseuse.


Techniquement, Fortran ne crée pas de tableau inutile. C'est le compilateur que vous utilisez qui crée le tableau temporaire. De plus, au lieu de deviner ce que fait gfortran, vous pouvez regarder les décharges intermédiaires. Essayez -fdump-tree-original et -fdump-tree-optimisé (le d de fin peut être incorrect). Enfin, au lieu d'essayer d'être intelligent, écrivez le code de la manière que vous désirez, c'est-à-dire r = 0, do i=1, size(xs); r = r + f(xs(i)); end do


@evets Peut-être vouliez-vous adresser cela à kdb?



2
votes

Pour expérimenter, j'ai essayé ce code (qui calcule la somme de l'inverse de arr ).

val = 0
do i = 1, num
    val = val + testfunc( arr(i) )
enddo

Ensuite, sur mon ordinateur (mac2012), "gfortran-10 -O2 test.f90 && time ./a.out" donne

 val =    18.997896413852555       in   0.225999996      (s)
 val =    18.997896413852555       in   0.252000004      (s)
 val =    18.997896413852555       in   0.246999994      (s)
 ...
real    0m3.909s
user    0m3.009s
sys 0m0.867s

et "gfortran-10 -O3" donne

 val =    18.997896413852555       in   0.437000006      (s)
 val =    18.997896413852555       in   0.453999996      (s)
 val =    18.997896413852555       in   0.437999994      (s)
 ...
real    0m5.946s
user    0m5.069s
sys 0m0.842s

Dans les deux cas, la commande htop affiche ~ 1,5 Go alloué, ce qui peut être raisonnable si les lignes L1 et L2 utilisent une matrice temporaire (chacune ~ 800 Mo avec ~ 0,3 s pour l'allocation).

Comme il n'y a pas de syntaxe pour créer des "itérateurs", j'ai essayé de faire de testfunc() un élémentaire (ou un élémentaire impur). La seule différence ici est les lignes marquées avec <-- .

program main
    use iso_fortran_env, only: dp => real64
    implicit none
    real(dp) val
    real(dp), allocatable :: arr(:)
    integer num, loop, i, t1, t2, trate

    num = 10**8

    arr = [( i, i = 1, num )]

    do loop = 1, 10
        call system_clock( t1 )

        val = sum( testfunc( arr ) )   !<--

        call system_clock( t2, trate )
        print *, "val = ", val, " in ", (t2 - t1) / real(trate), " (s)"
    enddo

contains

impure elemental &   !<--
function testfunc( x ) result( ret )
    real(dp), intent(in) :: x
    real(dp) :: ret

    ret = 1.0_dp / x
end

end program

Ensuite, "gfortran-10 -O2" donne

 val =    18.997896413852555       in   0.875000000      (s)
 val =    18.997896413852555       in   0.888000011      (s)
 val =    18.997896413852555       in   0.833000004      (s)
 ...

real    0m9.986s
user    0m5.738s
sys 0m4.210s

et "gfortran-10 -O3" donne

 val =    18.997896413852555       in    1.02999997      (s)
 val =    18.997896413852555       in    1.10099995      (s)
 val =    18.997896413852555       in    1.17600000      (s)
 ...

real    0m12.575s
user    0m8.142s
sys 0m4.387s

La commande htop affiche ~ 800 Mo, il semble donc que seul arr est alloué.

À titre de comparaison, le code suivant calcule val avec une boucle do explicite (en utilisant une version scalaire de testfunc() )

program main
    use iso_fortran_env, only: dp => real64
    implicit none
    real(dp) val
    real(dp), allocatable :: arr(:)
    integer num, loop, i, t1, t2, trate

    num = 10**8
    arr = [( i, i = 1, num )]   !! L1

    do loop = 1, 10
        call system_clock( t1 )

        val = sum( [( testfunc( arr(i) ), i = 1, num )] )   !! L2

        call system_clock( t2, trate )
        print *, "val = ", val, " in ", (t2 - t1) / real(trate), " (s)"
    enddo

contains

function testfunc( x ) result( ret )
    real(dp), intent(in) :: x
    real(dp) :: ret

    ret = 1.0_dp / x
end

end program

qui a donné le même timing avec le deuxième code avec elemental + testfunc(arr) ci-dessus (sur mon mac).

Dans tout le code ci-dessus, l'option -Ofast entraîné une erreur d'exécution ( Illegal instruction ). Mais cela était dû à la ligne L1 ( arr = [( i, i = 1, num )] ). Si -Ofast arr avance et que je le remplis avec une boucle explicite, -Ofast également fonctionné sans problème (donnant presque le même timing avec -O3 dans ce cas).


5 commentaires

Donc quelle est la réponse?


peut-être "il n'y a pas d'itérateur dans la norme" ...?


Alors pourquoi n'avez-vous pas mis ce i dans votre réponse? Peut-être que je suis juste en train de râler et que les gens votent pour votre réponse et ne votent pas pour la mienne, mais à mon avis, la réponse devrait d'abord et avant tout répondre à la question posée. Et puis montrez des expériences intéressantes sur les performances. Peut-être que je suis juste od et grincheux.


@VladimirF Salut, je pensais juste que ta (première) réponse donnait déjà la réponse rigoureuse (de la norme), ce qui, je pense, devrait être acceptée (cc: @kdd), j'ai donc ajouté quelques expériences sous un angle différent (car OP semblait s'intéresser aussi aux codes connexes). Quoi qu'il en soit, je pense que "élémentaire" est une fonctionnalité très utile.


OK, je pense qu'il est toujours bon de dire le point saillant d'une nouvelle réponse, même si cela a déjà été dit.



1
votes

Si vous créez la fonction calculateMomentum et élémentaire, elle peut être utilisée à la fois pour les valeurs scalaires et pour les tableaux

par exemple avec:

    type body
        real :: mass
    end type

vous pouvez appliquer ce qui précède à un tableau de v

integer::i
type(body) :: ball
real, allocatable :: v(:), p(:)
real :: tot_p
allocate(v(10))
v = [ (10+i, i=1, 10) ]
p = calculateMomentum(ball, v)
tot_p = sum(p)

avoir un tableau intermédiaire pour contenir les valeurs est avantageux car il maintient les données à proximité (probablement dans la ligne de cache) et la fonction sum() serait aussi rapide que possible.

imaginez le body du type comme suit par exemple

    elemental function calculateMomentum(obj, v) result(p)
    class(body), intent(in) :: obj
    real, intent(in) :: v
    real :: p
        p = obj%mass * v
    end function


0 commentaires