module opera1

  use adj_mod
  use tcor1
  implicit none
  integer, private				:: nw,nx,nz,nreal
  contains
!==================================================
  subroutine opera_init(nw_in,nx_in,nz_in,nreal_in)
    integer	:: nw_in,nx_in,nz_in,nreal_in
    nw=nw_in
    nx=nx_in
    nz=nz_in
    nreal=nreal_in
  end subroutine
!==================================================
  function opera_lop ( adj, add, xx, yy, x0) result(stat)
    integer            :: stat 
    logical,intent(in) :: adj,add 
    real,dimension(:)  :: xx,yy,x0
    call adjnull (adj,add,xx,yy )
    call opera_lop2(adj,add,xx,yy,x0 )
    stat=0
  end function 
!==================================================
  subroutine opera_lop2(adj,add,xx,yy,x0)
    logical,intent(in)   :: adj,add 
    logical              :: T=.true., F=.false.
    real, dimension (:)  :: xx 
    real, dimension (:)  :: yy 
    real, dimension (:)  :: x0 
    integer              :: i,j,k,stat
    integer, external    :: omp_get_thread_num,omp_get_num_threads

!!$OMP PARALLEL DO PRIVATE (inode,&
!                               j,&
!                            mini,&
!                            maxi,&
!                               i) schedule(dynamic,1)  num_threads(nodes) 



!!!! computar o yy da forma que eu quiser
    do j = 1, nx*nz
       inode = omp_get_thread_num()+1
       k = (j-1)*nx*nz
       do i = 1, nx*nz
          if ( adj) then
              stat = tcor1_lop ( T, F, xx((i-1)*nw+1:i*nw), yy, x0((j-1)*nw+1:j*nw))
          else
              stat = tcor1_lop ( F, F, xx((i-1)*nw+1:i*nw), yy, x0((j-1)*nw+1:j*nw))
!              if (i.eq.j) then
!                 yy=yy/maxval(abs(yy))
!              end if
          end if
       end do
    end do

!!$OMP END PARALLEL DO
  end subroutine 
end module opera1


