1! 2! 3! Fortran kernel for the MDot() vector routine 4! 5#include <petsc/finclude/petscsys.h> 6! 7pure subroutine FortranMDot4(x,y1,y2,y3,y4,n,sum1,sum2,sum3,sum4) 8 implicit none (type, external) 9 PetscScalar, intent(inout) :: sum1,sum2,sum3,sum4 10 PetscScalar, intent(in) :: x(*),y1(*),y2(*),y3(*),y4(*) 11 PetscInt, intent(in) :: n 12 13 PetscInt :: i 14 15 PETSC_AssertAlignx(16,x(1)) 16 PETSC_AssertAlignx(16,y1(1)) 17 PETSC_AssertAlignx(16,y2(1)) 18 PETSC_AssertAlignx(16,y3(1)) 19 PETSC_AssertAlignx(16,y4(1)) 20 21 do i=1,n 22 sum1 = sum1 + x(i)*PetscConj(y1(i)) 23 sum2 = sum2 + x(i)*PetscConj(y2(i)) 24 sum3 = sum3 + x(i)*PetscConj(y3(i)) 25 sum4 = sum4 + x(i)*PetscConj(y4(i)) 26 end do 27end subroutine FortranMDot4 28 29pure subroutine FortranMDot3(x,y1,y2,y3,n,sum1,sum2,sum3) 30 implicit none (type, external) 31 PetscScalar, intent(inout) :: sum1,sum2,sum3 32 PetscScalar, intent(in) :: x(*),y1(*),y2(*),y3(*) 33 PetscInt, intent(in) :: n 34 35 PetscInt :: i 36 37 PETSC_AssertAlignx(16,x(1)) 38 PETSC_AssertAlignx(16,y1(1)) 39 PETSC_AssertAlignx(16,y2(1)) 40 PETSC_AssertAlignx(16,y3(1)) 41 42 do i=1,n 43 sum1 = sum1 + x(i)*PetscConj(y1(i)) 44 sum2 = sum2 + x(i)*PetscConj(y2(i)) 45 sum3 = sum3 + x(i)*PetscConj(y3(i)) 46 end do 47end subroutine FortranMDot3 48 49pure subroutine FortranMDot2(x,y1,y2,n,sum1,sum2) 50 implicit none (type, external) 51 PetscScalar, intent(inout) :: sum1,sum2 52 PetscScalar, intent(in) :: x(*),y1(*),y2(*) 53 PetscInt, intent(in) :: n 54 55 PetscInt :: i 56 57 PETSC_AssertAlignx(16,x(1)) 58 PETSC_AssertAlignx(16,y1(1)) 59 PETSC_AssertAlignx(16,y2(1)) 60 61 do i=1,n 62 sum1 = sum1 + x(i)*PetscConj(y1(i)) 63 sum2 = sum2 + x(i)*PetscConj(y2(i)) 64 end do 65end subroutine FortranMDot2 66 67pure subroutine FortranMDot1(x,y1,n,sum1) 68 implicit none (type, external) 69 PetscScalar, intent(inout) :: sum1 70 PetscScalar, intent(in) :: x(*),y1(*) 71 PetscInt, intent(in) :: n 72 73 PetscInt :: i 74 75 PETSC_AssertAlignx(16,x(1)) 76 PETSC_AssertAlignx(16,y1(1)) 77 78 do i=1,n 79 sum1 = sum1 + x(i)*PetscConj(y1(i)) 80 end do 81 82end subroutine FortranMDot1 83