Actual source code: fmdot.F

  1: !
  2: !
  3: !    Fortran kernel for the MDot() vector routine
  4: !
 5:  #include include/finclude/petscdef.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

 13:       PetscInt i

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

 22:       return
 23:       end

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

 30:       PetscInt i

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

 38:       return
 39:       end

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

 46:       PetscInt i

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

 53:       return
 54:       end


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

 62:       PetscInt i

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

 68:       return
 69:       end