Как векторизовать эти циклы?

Question

Я имею дело с 4-кратным l oop. Это очень медленно из-за встроенной функции c exp () во внутренней l oop. Вот небольшой пример:

program main
  implicit none
  integer :: n, m, k, w
  real*8 :: a(100), b(100), c(1000), d(4)
  real*8 :: self(4,1000)

  a = 1.0d0
  b = 1.0d0
  c = 1.0d0
  d = 1.0d0

  self=0.0d0
  do n = 1, 100
    do m = 1, 100
      do k = 1, 1000
        do w = 1, 4
          self(w,k) = self(w,k) + exp( ( (c(k)-a(n))**2 + (c(k)-b(m))**2 ) / (2.0d0*d(w)**2) )
        enddo
      enddo
    enddo
  enddo

  ! my optimization:
  self=0.0d0
  do n = 1, 100
    do m = 1, 100
      !do k = 1, 1000
        do w = 1, 4
          self(w,:) = self(w,:) + exp( ( (c(:)-a(n))**2 + (c(:)-b(m))**2 ) / (2.0d0*d(w)**2) )
        enddo
      !enddo
    enddo
  enddo

end program

Похоже, что функция fortran intrinsi c exp () неэффективна. Тем не менее, я не хочу переписывать exp () как эквивалентные выражения.

Answer 1

Я разместил свой ответ здесь, основываясь на некоторых хороших комментариях выше:

program main
  implicit none
  integer :: n, m, w
  real*8 :: a(100), b(100), c(1000), d(4)
  real*8 :: self(1000,4), exp_tmp(1000,4,100)

  a = 1.0d0
  b = 1.0d0
  c = 1.0d0
  d = 1.0d0

  do n = 1, 100
    do w = 1, 4
      exp_tmp(:,w,n) = exp( 0.5d0 * ( (c(:)-a(n))/d(w) )**2 )
    enddo
  enddo

  self=0.0d0
  do n = 1, 100
    do m = 1, 100
      do w = 1, 4
        self(:,w) = self(:,w) + exp_tmp(:,w,n) * exp_tmp(:,w,m)
      enddo
    enddo
  enddo

end program

Я не знаю, можно ли векторизовать циклы, чтобы они были компактными.