Actual source code: fmaxpy_bgl.F

  1: ! ======= BGL Routines =====
  2: !
  3: !
  4: !    Fortran kernel for the MAXPY() vector routine
  5: !
 6:  #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:       call flush_(6)

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

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

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

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

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

 85:       return
 86:       end