# 1 "<stdin>"
# 1 "<built-in>"
# 1 "<command line>"
# 1 "<stdin>"
module allder
  use adj_mod
  implicit none
  integer, private :: nx
  integer, private :: ny
  real ,private, dimension (:), pointer :: pp
  contains
  subroutine allder_init ( nx_in,ny_in,pp_in )
    integer :: nx_in
    integer :: ny_in
    real , dimension (:), pointer :: pp_in
    nx = nx_in
    ny = ny_in
    pp => pp_in
  end subroutine
  function allder_lop ( adj, add, xx, yy) result(stat)
    integer :: stat
    logical,intent(in) :: adj,add
    real,dimension(:) :: xx,yy
    call adjnull (adj,add,xx,yy )
    call allder_lop2(adj,add,xx,yy )
    stat=0
  end function
  subroutine allder_lop2(adj,add,xx,yy)
    logical,intent(in) :: adj,add
    real, dimension (nx,ny) :: xx
    real, dimension (nx,ny) :: yy
    integer :: ix, iy,i2
    real :: a, ap, am, p
    integer, parameter :: dp = selected_real_kind(10,50)
    real(kind=dp) :: t
    do iy = 2,ny
      do ix = 2,nx-1
        p = pp(ix+(iy-2)*nx)
        ap = -3.-2.*p
        am = 3.-2.*p
        a = 4.*p
        ap=ap/12.
        am=am/12.
        a=a/12.
        if ( adj) then
          xx(ix-1,iy) = xx(ix-1,iy) + real(yy(ix,iy-1),dp) * real&
            &(am,dp)
          xx(ix-1,iy-1) = xx(ix-1,iy-1) - real(yy(ix,iy-1),dp) * real&
            &(ap,dp)
          xx(ix, iy-1) = xx(ix, iy-1) - real(yy(ix,iy-1),dp) * real&
            &(a,dp)
          xx(ix, iy) = xx(ix, iy) + yy(ix,iy) * a
          xx(ix+1,iy-1) = xx(ix+1,iy-1) - real(yy(ix,iy-1),dp) * real&
            &(am,dp)
          xx(ix+1,iy) = xx(ix+1,iy) +real( yy(ix,iy-1),dp) * real&
            &(ap,dp)
        else
          t= yy(ix,iy-1)
          t= t+real( (xx(ix-1,iy)-xx(ix+1,iy-1)),dp) *real( am,dp&
            &)
          t= t+real((xx(ix, iy)-xx(ix, iy-1)),dp) * real(a,dp)
          t= t+ real((xx(ix+1,iy)-xx(ix-1,iy-1)),dp) *real( ap &
            &,dp)
          yy(ix,iy-1)=t
        end if
      end do
    end do
  end subroutine
  subroutine allder_close()
  end subroutine
end module
