module triangle2_mod
! triangle smoothing in 2-D
  use adj_mod
  use triangle_mod
  implicit none
  integer, private :: b1 
  integer, private :: b2 
  integer, private :: n1 
  integer, private :: n2 
  real ,private,  dimension (:,:), allocatable  :: tmp 
  contains
  subroutine triangle2_init ( b1_in,b2_in,n1_in,n2_in )
    integer    :: b1_in 
    integer    :: b2_in 
    integer    :: n1_in 
    integer    :: n2_in 
    b1 = b1_in 
    b2 = b2_in 
    n1 = n1_in 
    n2 = n2_in 
    if (.not. allocated(tmp)) then
      allocate(tmp ( n1,n2)) 
    end if
  end subroutine
  function triangle2_lop ( adj, add, x, y) result(stat)
    integer            :: stat 
    logical,intent(in) :: adj,add 
    real,dimension(:)  :: x,y 
    call adjnull (adj,add,x,y )
    call triangle2_lop2(adj,add,x,y )
    stat=0
  end function 
  subroutine triangle2_lop2(adj,add,x,y)
    logical,intent(in) :: adj,add 
    real, dimension (n1,n2)  :: x 
    real, dimension (n1,n2)  :: y 
    integer i1, i2, stat1, stat2
    if (b1 > 1) then
      call triangle_init (b1,n1)
      do i2 = 1, n2  
        if (b2 > 1) then
          if (adj) then
            stat1 = triangle_lop (adj, .false., tmp (:,i2), y (:,i2))
          else
            stat1 = triangle_lop (adj, .false., x (:,i2), tmp (:,i2))
          end if
        else
          stat1 = triangle_lop (adj, .true., x (:,i2), y (:,i2))
        end if
      end do 
      call triangle_close ()
    end if
    if (b2 > 1) then
      call triangle_init (b2,n2)
      do i1 = 1, n1  
        if (b1 > 1) then
          if (adj) then
            stat2 = triangle_lop (adj, .true., x (i1,:), tmp (i1,:))
          else
            stat2 = triangle_lop (adj, .true., tmp (i1,:), y (i1,:))
          end if
        else
          stat1 = triangle_lop (adj, .true., x (i1,:), y (i1,:))
        end if
      end do 
      call triangle_close ()
    end if
  end subroutine 
  subroutine triangle2_close()
    deallocate( tmp )
  end subroutine 
end module 

