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 end 28 29 subroutine FortranMDot3(x,y1,y2,y3,n,sum1,sum2,sum3) 30 implicit none 31 PetscScalar sum1,sum2,sum3 32 PetscScalar x(*),y1(*),y2(*),y3(*) 33 PetscInt n 34 PetscInt i 35 36 PETSC_AssertAlignx(16,x(1)) 37 PETSC_AssertAlignx(16,y1(1)) 38 PETSC_AssertAlignx(16,y2(1)) 39 PETSC_AssertAlignx(16,y3(1)) 40 do 10,i=1,n 41 sum1 = sum1 + x(i)*PetscConj(y1(i)) 42 sum2 = sum2 + x(i)*PetscConj(y2(i)) 43 sum3 = sum3 + x(i)*PetscConj(y3(i)) 44 10 continue 45 46 end 47 48 subroutine FortranMDot2(x,y1,y2,n,sum1,sum2) 49 implicit none 50 PetscScalar sum1,sum2,x(*),y1(*),y2(*) 51 PetscInt n 52 PetscInt i 53 54 PETSC_AssertAlignx(16,x(1)) 55 PETSC_AssertAlignx(16,y1(1)) 56 PETSC_AssertAlignx(16,y2(1)) 57 do 10,i=1,n 58 sum1 = sum1 + x(i)*PetscConj(y1(i)) 59 sum2 = sum2 + x(i)*PetscConj(y2(i)) 60 10 continue 61 62 end 63 64 subroutine FortranMDot1(x,y1,n,sum1) 65 implicit none 66 PetscScalar sum1,x(*),y1(*) 67 PetscInt n 68 PetscInt i 69 70 PETSC_AssertAlignx(16,x(1)) 71 PETSC_AssertAlignx(16,y1(1)) 72 do 10,i=1,n 73 sum1 = sum1 + x(i)*PetscConj(y1(i)) 74 10 continue 75 76 end 77