Actual source code: fmdot.F

  1: !
  2: !  "$Id: fmdot.F,v 1.10 2001/08/07 03:05:24 balay Exp $";
  3: !
  4: !    Fortran kernel for the MDot() vector routine
  5: !
 6:  #include include/finclude/petscdef.h
  7: !
  8:       subroutine FortranMDot4(x,y1,y2,y3,y4,n,sum1,sum2,sum3,sum4)
  9:       implicit none
 10:       PetscScalar  sum1,sum2,sum3,sum4
 11:       PetscScalar  x(*),y1(*),y2(*),y3(*),y4(*)
 12:       integer n

 14:       integer i

 16:       do 10,i=1,n
 17:         sum1 = sum1 + x(i)*PetscConj(y1(i))
 18:         sum2 = sum2 + x(i)*PetscConj(y2(i))
 19:         sum3 = sum3 + x(i)*PetscConj(y3(i))
 20:         sum4 = sum4 + x(i)*PetscConj(y4(i))
 21:  10   continue

 23:       return
 24:       end

 26:       subroutine FortranMDot3(x,y1,y2,y3,n,sum1,sum2,sum3)
 27:       implicit none
 28:       PetscScalar  sum1,sum2,sum3,x(*),y1(*),y2(*),y3(*)
 29:       integer n

 31:       integer i

 33:       do 10,i=1,n
 34:         sum1 = sum1 + x(i)*PetscConj(y1(i))
 35:         sum2 = sum2 + x(i)*PetscConj(y2(i))
 36:         sum3 = sum3 + x(i)*PetscConj(y3(i))
 37:  10   continue

 39:       return
 40:       end

 42:       subroutine FortranMDot2(x,y1,y2,n,sum1,sum2)
 43:       implicit none
 44:       PetscScalar  sum1,sum2,x(*),y1(*),y2(*)
 45:       integer n

 47:       integer i

 49:       do 10,i=1,n
 50:         sum1 = sum1 + x(i)*PetscConj(y1(i))
 51:         sum2 = sum2 + x(i)*PetscConj(y2(i))
 52:  10   continue

 54:       return
 55:       end


 58:       subroutine FortranMDot1(x,y1,n,sum1)
 59:       implicit none
 60:       PetscScalar  sum1,x(*),y1(*)
 61:       integer n

 63:       integer i

 65:       do 10,i=1,n
 66:         sum1 = sum1 + x(i)*PetscConj(y1(i))
 67:  10   continue

 69:       return
 70:       end