the power function is serialized under Intel Fortran

[mehdi__chinoune]

I tried to use Intel Advisor to vectorise my program. and it tells me that the compiler con not vectorise some loops because it contain _powr8i4 which correspond to x**i (x:double precision and i:integer).

Why did not intel Fortran Compiler come with a vectorized version of _powr8i4 (or generally power function)?

and How to vectorise it (other than using MKL)?

NB: I did use the most aggressive optimization flags : /O3 /QxHost /Qipo /fp:fast=2

[jimdempseyatthecove]

For + powers, try:

    module foo
    contains
    !dir$ attributes vector : vectorlength(4) :: vpowr8i4
    elemental function vpowr8i4(baseIN, expIN) result(ret)
        implicit none
        real(8), intent(in) :: baseIN
        integer(4), intent(in) :: expIN
        real(8) :: base, ret
        integer(4) :: exp
        base = baseIn
        exp = expIN
        ret = 1.0_8
        do while(exp .ne. 0)
            if(IAND(exp,1) .ne. 0) ret = ret * base;
            exp = ishft(exp, -1)
            base = base * base
        end do
    end function vpowr8i4
    end module foo
    
    program test_vpowr8i4
    use foo
    implicit none
    integer, parameter :: arraySize = 1024
    !dir$ attributes align : 64 :: array,arrayOut
    real(8) :: array(arraySize)
    real(8) :: arrayOut(arraySize)
    integer :: i,exp
    do i=1,arraySize
        array(i) = i
    end do
    print *,"Enter exponent"
    read *,exp
    print *,"Break on this line in Release build"
    !dir$ simd
    do i=1,arraySize
        arrayOut(i) = vpowr8i4(array(i), exp)
    end do
    
    do i=1,10
        print *, array(i), arrayOut(i)
    end do
    
    
    end program test_vpowr8i4

 Enter exponent
3
 Break on this line in Release build
   1.00000000000000        1.00000000000000
   2.00000000000000        8.00000000000000
   3.00000000000000        27.0000000000000
   4.00000000000000        64.0000000000000
   5.00000000000000        125.000000000000
   6.00000000000000        216.000000000000
   7.00000000000000        343.000000000000
   8.00000000000000        512.000000000000
   9.00000000000000        729.000000000000
   10.0000000000000        1000.00000000000
Press any key to continue . . .

Jim Dempsey

[mehdi__chinoune]

jimdempseyatthecove wrote:

For + powers, try:

    module foo
    contains
    !dir$ attributes vector : vectorlength(4) :: vpowr8i4
    elemental function vpowr8i4(baseIN, expIN) result(ret)
        implicit none
        real(8), intent(in) :: baseIN
        integer(4), intent(in) :: expIN
        real(8) :: base, ret
        integer(4) :: exp
        base = baseIn
        exp = expIN
        ret = 1.0_8
        do while(exp .ne. 0)
            if(IAND(exp,1) .ne. 0) ret = ret * base;
            exp = ishft(exp, -1)
            base = base * base
        end do
    end function vpowr8i4
    end module foo
    
    program test_vpowr8i4
    use foo
    implicit none
    integer, parameter :: arraySize = 1024
    !dir$ attributes align : 64 :: array,arrayOut
    real(8) :: array(arraySize)
    real(8) :: arrayOut(arraySize)
    integer :: i,exp
    do i=1,arraySize
        array(i) = i
    end do
    print *,"Enter exponent"
    read *,exp
    print *,"Break on this line in Release build"
    !dir$ simd
    do i=1,arraySize
        arrayOut(i) = vpowr8i4(array(i), exp)
    end do
    
    do i=1,10
        print *, array(i), arrayOut(i)
    end do
    
    
    end program test_vpowr8i4

 Enter exponent
3
 Break on this line in Release build
   1.00000000000000        1.00000000000000
   2.00000000000000        8.00000000000000
   3.00000000000000        27.0000000000000
   4.00000000000000        64.0000000000000
   5.00000000000000        125.000000000000
   6.00000000000000        216.000000000000
   7.00000000000000        343.000000000000
   8.00000000000000        512.000000000000
   9.00000000000000        729.000000000000
   10.0000000000000        1000.00000000000
Press any key to continue . . .

Jim Dempsey

Thank you, I know these tricks but I don't want to modify my code or to use extensions. I just asked if there is something i missed like flags or link ...