Actual source code: fmaxpy_bgl.F

  1: ! ======= BGL Routines =====
  2: !
  3: !
  4: !    Fortran kernel for the MAXPY() vector routine
  5: !
 6:  #include include/finclude/petscdef.h
  7: !

  9:       Subroutine FortranMAXPY4_BGL(x, a0, a1, a2, a3, y0, y1, y2, y3, n)
 10:       implicit none
 11:       PetscScalar  a0,a1,a2,a3
 12:       PetscScalar  x(*),y0(*)
 13:       PetscScalar y1(*),y2(*),y3(*)
 14:       PetscInt n, i

 16:       write(6,*) 'n=', n
 17:       call flush_(6)

 19:       call ALIGNX(16,x(1))
 20:       call ALIGNX(16,y0(1))
 21:       call ALIGNX(16,y1(1))
 22:       call ALIGNX(16,Y2(1))
 23:       call ALIGNX(16,Y3(1))
 24:       do i=1,n
 25:           x(i)  = x(i) + a0*y0(i) + a1*y1(i) + a2*y2(i) + a3*y3(i)
 26:       enddo
 27:       return
 28:       end

 30:       subroutine FortranMAXPY3_BGL(x,a0,a1,a2,y0,y1,y2,n)
 31:       implicit none
 32:       PetscScalar  a0,a1,a2,x(*)
 33:       PetscScalar y0(*),y1(*),y2(*)
 34:       PetscInt n
 35:       PetscInt i
 36:       call ALIGNX(16,x(1))
 37:       call ALIGNX(16,y0(1))
 38:       call ALIGNX(16,y1(1))
 39:       call ALIGNX(16,y2(1))
 40:       do 10,i=1,n
 41:          x(i) = x(i) + a0*y0(i) + a1*y1(i) + a2*y2(i)
 42:   10  continue
 43:       return
 44:       end

 46:       Subroutine FortranMAXPY2_BGL(x, a0, a1, y0, y1, n)
 47:       implicit none
 48:       PetscScalar  a0,a1,x(*)
 49:       PetscScalar  y0(*),y1(*)
 50:       PetscInt n, i
 51:       call ALIGNX(16,x(1))
 52:       call ALIGNX(16,y0(1))
 53:       call ALIGNX(16,y1(1))
 54:       do i=1,n
 55:           x(i)  = x(i) + a0*y0(i) + a1*y1(i)
 56:       enddo
 57:       return
 58:       end

 60:       subroutine Fortranxtimesy_BGL(x,y,z,n)
 61:       implicit none
 62:       PetscScalar  x(*),y(*),z(*)
 63:       PetscInt n
 64:       PetscInt i
 65:       call ALIGNX(16,x(1))
 66:       call ALIGNX(16,y(1))
 67:       call ALIGNX(16,z(1))
 68:       do 10,i=1,n
 69:          z(i) = x(i) * y(i)
 70:   10  continue
 71:       return
 72:       end

 74:       subroutine FortranAYPX_BGL(n,a,x,y)
 75:       implicit none
 76:       PetscScalar  a
 77:       PetscScalar  x(*),y(*)
 78:       PetscInt n
 79:       PetscInt i
 80:       call ALIGNX(16,x(1))
 81:       call ALIGNX(16,y(1))
 82:       do 10,i=1,n
 83:         y(i) = x(i) + a*y(i)
 84:  10   continue

 86:       return
 87:       end