openlang

простые циклы очень медленные, gfortran, массив 200000 x 64

Небольшая утилита для калибровки, код ниже. Результат правильный, но слишком медленно. detrend+clip занимает несколько секунд, rms тоже медленнее, чем должен быть, при том что на элемент всего одно умножение и одно сложение. Пробовал -O3 и -march=native, почти ничего не изменилось. Подозревал sum и real(dp) на Apple silicon. Эти две медленные стадии — одна проблема или две? Не понимаю, почему это столько стоит.

Makefile

FC = gfortran
FFLAGS = -O3 -march=native -g

sensor_cal: src/sensor_mod.f90 src/main.f90
	$(FC) $(FFLAGS) -o $@ $^

clean:
	rm -f sensor_cal *.mod

src/sensor_mod.f90

module sensor_mod
  implicit none
  integer, parameter :: dp = kind(1.0d0)
contains

  ! remove the mean of one record (all channels at one time step)
  subroutine detrend_record(v, n)
    integer, intent(in) :: n
    real(dp), intent(inout) :: v(n)
    real(dp) :: m
    m = sum(v) / n
    v = v - m
  end subroutine detrend_record

  subroutine clip_record(v, n, lim)
    integer, intent(in) :: n
    real(dp), intent(inout) :: v(n)
    real(dp), intent(in) :: lim
    integer :: k
    do k = 1, n
      if (v(k) > lim) v(k) = lim
      if (v(k) < -lim) v(k) = -lim
    end do
  end subroutine clip_record

  function channel_rms(a, nrec, nchan) result(r)
    integer, intent(in) :: nrec, nchan
    real(dp), intent(in) :: a(nrec, nchan)
    real(dp) :: r(nchan)
    integer :: i, j
    r = 0.0_dp
    do i = 1, nrec
      do j = 1, nchan
        r(j) = r(j) + a(i, j)**2
      end do
    end do
    r = sqrt(r / nrec)
  end function channel_rms

end module sensor_mod

src/main.f90

program sensor_cal
  use sensor_mod
  implicit none
  integer :: nrec, nchan, i, t0, t1, rate
  real(dp), allocatable :: readings(:,:), rms(:)
  character(len=32) :: arg
  logical :: verbose

  verbose = .false.
  nrec = 200000
  nchan = 64

  if (command_argument_count() >= 1) then
    call get_command_argument(1, arg)
    read(arg, *) nrec
  end if
  if (command_argument_count() >= 2) then
    call get_command_argument(2, arg)
    read(arg, *) nchan
  end if
  if (command_argument_count() >= 3) verbose = .true.

  allocate(readings(nrec, nchan), rms(nchan))
  call random_number(readings)
  readings = 100.0_dp * readings

  call system_clock(t0, rate)
  do i = 1, nrec
    call detrend_record(readings(i, :), nchan)
    call clip_record(readings(i, :), nchan, 40.0_dp)
    if (verbose .and. mod(i, 50000) == 0) print *, 'record', i
  end do
  call system_clock(t1)
  print '(a,f8.3,a)', 'detrend+clip: ', real(t1 - t0) / rate, ' s'

  call system_clock(t0)
  rms = channel_rms(readings, nrec, nchan)
  call system_clock(t1)
  print '(a,f8.3,a)', 'rms:          ', real(t1 - t0) / rate, ' s'

  print '(a,3f10.4)', 'first rms values: ', rms(1:3)

  ! TODO: write rms to calibration file, format not agreed with the lab yet
end program sensor_cal

Ассемблер я не смотрел и предупреждения не включал. Тайминги выше по памяти, не точные.

Переведено с английскогоПоказать оригинал

Разговор продолжается 2

Следуя доводам
  1. Ссылка на ответ

    Одна проблема, не две. Fortran хранит массивы по столбцам (column-major), поэтому readings(i,:) идёт с шагом и при каждом вызове копируется во временный массив. В rms поменяйте циклы местами. Лучше объявите массив как (nchan, nrec). Скомпилируйте с -Warray-temporaries, чтобы увидеть копии.

    Переведено с английскогоПоказать оригинал

Другие стороны того же вопроса

  1. running_mean тормозит, время уходит на memmove, а не на арифметику

    Сглаживаю данные датчиков. Массив readings(nstation, ntime), станция первым индексом, потому что так устроен файл с данными и…

    Ответы 1Последняя активность fortran
  2. Какое время сборки вы приводите?

    Каждая тема про разбиение крейтов заканчивается тем, что кто-то выкладывает одно время сборки. Чисел два, и меняются они в ра…

    Ответы 2Последняя активность rust
  3. Нужен ли этой очереди сохранённый счётчик?

    Вызывающий код спрашивает у этой очереди только List.length items > 8, зато спрашивает постоянно. Очередь представлена одним…

    Ответы 2Последняя активность ocaml