matvec_struct.f File Reference

Go to the source code of this file.

Functions/Subroutines
subroutine	matvec_struct (n, x, y, nelt, ia, ja, a, isym)

Function/Subroutine Documentation

matvec_struct()

subroutine matvec_struct	(	integer	n,
		real*8, dimension(*)	x,
		real*8, dimension(*)	y,
		integer	nelt,
		integer, dimension(*)	ia,
		integer, dimension(*)	ja,
		real*8, dimension(*)	a,
		integer	isym
	)

      use omp_lib
!
      implicit none
!
      integer ia(*),ja(*),i,j,l,n,nelt,isym,nd,na
      real*8 y(*),x(*),a(*)
!
!     number of off-diagonal terms
!
      nd=nelt-n
!
      if(isym.eq.0) then
         na=nd/2
!
!        non-symmetric
!
!        diagonal terms
!
c$omp parallel default(none)
c$omp& shared(n,x,a,y,nd,ja,ia,na)
c$omp& private(i,l,j)
c$omp do
         do i=1,n
            y(i)=a(nd+i)*x(i)
         enddo
c$omp end do
!     
!        off-diagonal terms
!     
!        number of upper triangular terms
!
c$omp do
         do j=1,n
            do l=ja(j),ja(j+1)-1
               i=ia(l)
c$omp atomic
               y(i)=y(i)+a(l)*x(j)
            enddo
         enddo
c$omp end do
!
c$omp do
         do j=1,n
            do l=ja(j),ja(j+1)-1
               i=ia(l)
               y(j)=y(j)+a(l+na)*x(i)
            enddo
         enddo
c$omp end do
c$omp end parallel
!
      else
!
!        symmetric
!     
!        diagonal terms
!
c$omp parallel default(none)
c$omp& shared(n,x,a,y,nd,ja,ia)
c$omp& private(i,l,j)
c$omp do
         do i=1,n
            y(i)=a(nd+i)*x(i)
         enddo
c$omp end do
!
!        off-diagonal terms
!     
c$omp do
         do j=1,n
            do l=ja(j),ja(j+1)-1
               i=ia(l)
c$omp atomic
               y(i)=y(i)+a(l)*x(j)
            enddo
         enddo
c$omp end do
!
c$omp do
         do j=1,n
            do l=ja(j),ja(j+1)-1
               i=ia(l)
               y(j)=y(j)+a(l)*x(i)
            enddo
         enddo
c$omp end do
c$omp end parallel
      endif
!
      return