1c96caaccSSatish Balay! 2c96caaccSSatish Balay! 3c96caaccSSatish Balay! Fortran kernel for the MDot() vector routine 4c96caaccSSatish Balay! 5c96caaccSSatish Balay#include <petsc/finclude/petscsys.h> 6c96caaccSSatish Balay! 70ccf82acSMartin Diehlpure subroutine FortranMDot4(x,y1,y2,y3,y4,n,sum1,sum2,sum3,sum4) 80ccf82acSMartin Diehl implicit none (type, external) 90ccf82acSMartin Diehl PetscScalar, intent(inout) :: sum1,sum2,sum3,sum4 100ccf82acSMartin Diehl PetscScalar, intent(in) :: x(*),y1(*),y2(*),y3(*),y4(*) 110ccf82acSMartin Diehl PetscInt, intent(in) :: n 120113e719SMartin Diehl 13*d66e387eSMartin Diehl PetscInt :: i 14c96caaccSSatish Balay 15c96caaccSSatish Balay PETSC_AssertAlignx(16,x(1)) 16c96caaccSSatish Balay PETSC_AssertAlignx(16,y1(1)) 17c96caaccSSatish Balay PETSC_AssertAlignx(16,y2(1)) 18c96caaccSSatish Balay PETSC_AssertAlignx(16,y3(1)) 19c96caaccSSatish Balay PETSC_AssertAlignx(16,y4(1)) 20c96caaccSSatish Balay 210113e719SMartin Diehl do i=1,n 22c96caaccSSatish Balay sum1 = sum1 + x(i)*PetscConj(y1(i)) 23c96caaccSSatish Balay sum2 = sum2 + x(i)*PetscConj(y2(i)) 24c96caaccSSatish Balay sum3 = sum3 + x(i)*PetscConj(y3(i)) 25c96caaccSSatish Balay sum4 = sum4 + x(i)*PetscConj(y4(i)) 260113e719SMartin Diehl end do 270113e719SMartin Diehlend subroutine FortranMDot4 28c96caaccSSatish Balay 290ccf82acSMartin Diehlpure subroutine FortranMDot3(x,y1,y2,y3,n,sum1,sum2,sum3) 300ccf82acSMartin Diehl implicit none (type, external) 310ccf82acSMartin Diehl PetscScalar, intent(inout) :: sum1,sum2,sum3 320ccf82acSMartin Diehl PetscScalar, intent(in) :: x(*),y1(*),y2(*),y3(*) 330ccf82acSMartin Diehl PetscInt, intent(in) :: n 340113e719SMartin Diehl 35*d66e387eSMartin Diehl PetscInt :: i 36c96caaccSSatish Balay 37c96caaccSSatish Balay PETSC_AssertAlignx(16,x(1)) 38c96caaccSSatish Balay PETSC_AssertAlignx(16,y1(1)) 39c96caaccSSatish Balay PETSC_AssertAlignx(16,y2(1)) 40c96caaccSSatish Balay PETSC_AssertAlignx(16,y3(1)) 410113e719SMartin Diehl 420113e719SMartin Diehl do i=1,n 43c96caaccSSatish Balay sum1 = sum1 + x(i)*PetscConj(y1(i)) 44c96caaccSSatish Balay sum2 = sum2 + x(i)*PetscConj(y2(i)) 45c96caaccSSatish Balay sum3 = sum3 + x(i)*PetscConj(y3(i)) 460113e719SMartin Diehl end do 470113e719SMartin Diehlend subroutine FortranMDot3 48c96caaccSSatish Balay 490ccf82acSMartin Diehlpure subroutine FortranMDot2(x,y1,y2,n,sum1,sum2) 500ccf82acSMartin Diehl implicit none (type, external) 510ccf82acSMartin Diehl PetscScalar, intent(inout) :: sum1,sum2 520ccf82acSMartin Diehl PetscScalar, intent(in) :: x(*),y1(*),y2(*) 530ccf82acSMartin Diehl PetscInt, intent(in) :: n 540113e719SMartin Diehl 55*d66e387eSMartin Diehl PetscInt :: i 56c96caaccSSatish Balay 57c96caaccSSatish Balay PETSC_AssertAlignx(16,x(1)) 58c96caaccSSatish Balay PETSC_AssertAlignx(16,y1(1)) 59c96caaccSSatish Balay PETSC_AssertAlignx(16,y2(1)) 600113e719SMartin Diehl 610113e719SMartin Diehl do i=1,n 62c96caaccSSatish Balay sum1 = sum1 + x(i)*PetscConj(y1(i)) 63c96caaccSSatish Balay sum2 = sum2 + x(i)*PetscConj(y2(i)) 640113e719SMartin Diehl end do 650113e719SMartin Diehlend subroutine FortranMDot2 66c96caaccSSatish Balay 670ccf82acSMartin Diehlpure subroutine FortranMDot1(x,y1,n,sum1) 680ccf82acSMartin Diehl implicit none (type, external) 690ccf82acSMartin Diehl PetscScalar, intent(inout) :: sum1 700ccf82acSMartin Diehl PetscScalar, intent(in) :: x(*),y1(*) 710ccf82acSMartin Diehl PetscInt, intent(in) :: n 720113e719SMartin Diehl 73*d66e387eSMartin Diehl PetscInt :: i 74c96caaccSSatish Balay 75c96caaccSSatish Balay PETSC_AssertAlignx(16,x(1)) 76c96caaccSSatish Balay PETSC_AssertAlignx(16,y1(1)) 77c96caaccSSatish Balay 780113e719SMartin Diehl do i=1,n 790113e719SMartin Diehl sum1 = sum1 + x(i)*PetscConj(y1(i)) 800113e719SMartin Diehl end do 810113e719SMartin Diehl 820113e719SMartin Diehlend subroutine FortranMDot1 83