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