! Horizontal basic numerical operations for 2-dimensional arrays.
! These are the basic building blocks for horizontal compact and
! conventional differencing, quadrature, midpoint interpolation
! and filtering.
!############# 2-dimensions, index-i active:
subroutine ud1c2i(c,d,bnd1, & 1
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: c
real, dimension(ims:ime , jms:jme), &
INTENT(OUT ) :: d
real, INTENT(IN ) :: bnd1
integer :: i,j
do j=jts,jte
do i=its,ite
d(i,j)=bnd1*(c(i+1,j)-c(i-1,j))
enddo
enddo
end subroutine ud1c2i
subroutine ud2c2i(c,d,bnd1,bnd2, & 1
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: c
real, dimension(ims:ime , jms:jme), &
INTENT(OUT ) :: d
real, INTENT(IN ) :: bnd1,bnd2
integer :: i,j
do j=jts,jte
do i=its,ite
d(i,j)=bnd1*(c(i+1,j)-c(i-1,j))+bnd2*(c(i+2,j)-c(i-2,j))
enddo
enddo
end subroutine ud2c2i
subroutine sa1a2i(a,bnm1, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: bnm1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=bnm1*(a(i,j)+a(i+1,j))
enddo
enddo
end subroutine sa1a2i
subroutine sq1a2i(a,bnq, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: bnq
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i-1,j)+bnq*a(i,j)
enddo
enddo
end subroutine sq1a2i
subroutine sd1b2i(a,bnqi, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: bnqi
integer :: i,j
do j=jts,jte
do i=ite,its,-1
a(i,j)=bnqi*(a(i,j)-a(i-1,j))
enddo
enddo
end subroutine sd1b2i
subroutine ef1a2i(a,knq1, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: knq1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)+knq1*a(i+1,j)
enddo
enddo
end subroutine ef1a2i
subroutine ef2a2i(a,knq1,knq2, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: knq1,knq2
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)+knq1*a(i+1,j)+knq2*a(i+2,j)
enddo
enddo
end subroutine ef2a2i
subroutine ef1b2i(a,knq1, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: knq1
integer :: i,j
do j=jts,jte
do i=ite,its,-1
a(i,j)=a(i,j)+knq1*a(i-1,j)
enddo
enddo
end subroutine ef1b2i
subroutine ef2b2i(a,knq1,knq2, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: knq1,knq2
integer :: i,j
do j=jts,jte
do i=ite,its,-1
a(i,j)=a(i,j)+knq1*a(i-1,j)+knq2*a(i-2,j)
enddo
enddo
end subroutine ef2b2i
subroutine rf1a2i(a,lnq1, & 8
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i-1,j)
enddo
enddo
end subroutine rf1a2i
subroutine rf2a2i(a,lnq1,lnq2, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i-1,j)-lnq2*a(i-2,j)
enddo
enddo
end subroutine rf2a2i
subroutine rf1b2i(a,lnq1, & 8
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1
integer :: i,j
do j=jts,jte
do i=ite,its,-1
a(i,j)=a(i,j)-lnq1*a(i+1,j)
enddo
enddo
end subroutine rf1b2i
subroutine rf2b2i(a,lnq1,lnq2, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2
integer :: i,j
do j=jts,jte
do i=ite,its,-1
a(i,j)=a(i,j)-lnq1*a(i+1,j)-lnq2*a(i+2,j)
enddo
enddo
end subroutine rf2b2i
subroutine ra1a2i(a,p,lnq1,fnq, & 8
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(1, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In rala2i; fnq exceeds tile width'
do j=jts,jte
p(1,j)=0.
enddo
do j=jts,jte
do i=ite+1-fnq,ite
p(1,j)=a(i,j)-lnq1*p(1,j)
enddo
enddo
end subroutine ra1a2i
subroutine ra2a2i(a,p,lnq1,lnq2,fnq, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(2, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In ra2a2i; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ra2a2i; fnq is not even'
do j=jts,jte
do i=1,2
p(i,j)=0.
enddo
enddo
do j=jts,jte
do i=ite+1-fnq,ite,2
p(2,j)=a(i ,j)-lnq1*p(1,j)-lnq2*p(2,j)
p(1,j)=a(i+1,j)-lnq1*p(2,j)-lnq2*p(1,j)
enddo
enddo
end subroutine ra2a2i
subroutine ra1b2i(a,p,lnq1,fnq, & 8
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(1, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In ralb2i; fnq exceeds tile width'
do j=jts,jte
p(1,j)=0.
enddo
do j=jts,jte
do i=its-1+fnq,its,-1
p(1,j)=a(i,j)-lnq1*p(1,j)
enddo
enddo
end subroutine ra1b2i
subroutine ra2b2i(a,p,lnq1,lnq2,fnq, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(2, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In ra2b2i; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ra2b2i; fnq is not even'
do j=jts,jte
do i=1,2
p(i,j)=0.
enddo
enddo
do j=jts,jte
do i=its-1+fnq,its,-2
p(2,j)=a(i ,j)-lnq1*p(1,j)-lnq2*p(2,j)
p(1,j)=a(i-1,j)-lnq1*p(2,j)-lnq2*p(1,j)
enddo
enddo
end subroutine ra2b2i
subroutine mf11a2i(a,lnq1,knq1, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,knq1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i-1,j)+knq1*a(i+1,j)
enddo
enddo
end subroutine mf11a2i
subroutine mf12a2i(a,lnq1,knq1,knq2, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,knq1,knq2
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i-1,j) &
+knq1*a(i+1,j)+knq2*a(i+2,j)
enddo
enddo
end subroutine mf12a2i
subroutine mf21a2i(a,lnq1,lnq2,knq1, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2,knq1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i-1,j)-lnq2*a(i-2,j)+knq1*a(i+1,j)
enddo
enddo
end subroutine mf21a2i
subroutine mf22a2i(a,lnq1,lnq2,knq1,knq2, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2,knq1,knq2
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i-1,j)-lnq2*a(i-2,j) &
+knq1*a(i+1,j)+knq2*a(i+2,j)
enddo
enddo
end subroutine mf22a2i
subroutine ma11a2i(a,p,lnq1,knq1,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(1, jms:jme),&
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,knq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In mal1a2i; fnq exceeds tile width'
do j=jts,jte
p(1,j)=0.
enddo
do j=jts,jte
do i=ite+1-fnq,ite
p(1,j)=a(i,j)-lnq1*p(1,j)+knq1*a(i+1,j)
enddo
enddo
end subroutine ma11a2i
subroutine ma12a2i(a,p,lnq1,knq1,knq2,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(1, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,knq1,knq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In ma12a2i; fnq exceeds tile width'
do j=jts,jte
p(1,j)=0.
enddo
do j=jts,jte
do i=ite+1-fnq,ite
p(1,j)=a(i,j)-lnq1*p(1, j) &
+knq1*a(i+1,j)+knq2*a(i+2,j)
enddo
enddo
end subroutine ma12a2i
subroutine ma21a2i(a,p,lnq1,lnq2,knq1,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(2, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2,knq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In ma21a2i; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ma21a2i; fnq is not even'
do j=jts,jte
p(1,j)=0.
p(2,j)=0.
enddo
do j=jts,jte
do i=ite+1-fnq,ite,2
p(2,j)=a(i, j)-lnq1*p(1, j)-lnq2*p(2,j) &
+knq1*a(i+1,j)
p(1,j)=a(i+1,j)-lnq1*p(2, j)-lnq2*p(1,j) &
+knq1*a(i+2,j)
enddo
enddo
end subroutine ma21a2i
subroutine ma22a2i(a,p,lnq1,lnq2,knq1,knq2,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(2, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2,knq1,knq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In ma22a2i; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ma22a2i; fnq is not even'
do j=jts,jte
p(1,j)=0.
p(2,j)=0.
enddo
do j=jts,jte
do i=ite+1-fnq,ite,2
p(2,j)=a(i, j)-lnq1*p(1, j)-lnq2*p(2, j) &
+knq1*a(i+1,j)+knq2*a(i+2,j)
p(1,j)=a(i+1,j)-lnq1*p(2, j)-lnq2*p(1, j) &
+knq1*a(i+2,j)+knq2*a(i+3,j)
enddo
enddo
end subroutine ma22a2i
subroutine mf11b2i(a,lnq1,knq1, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,knq1
integer :: i,j
do j=jts,jte
do i=ite,its,-1
a(i,j)=a(i,j)-lnq1*a(i+1,j)+knq1*a(i-1,j)
enddo
enddo
end subroutine mf11b2i
subroutine mf12b2i(a,lnq1,knq1,knq2, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,knq1,knq2
integer :: i,j
do j=jts,jte
do i=ite,its,-1
a(i,j)=a(i,j)-lnq1*a(i+1,j) &
+knq1*a(i-1,j)+knq2*a(i-2,j)
enddo
enddo
end subroutine mf12b2i
subroutine mf21b2i(a,lnq1,lnq2,knq1, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2,knq1
integer :: i,j
do j=jts,jte
do i=ite,its,-1
a(i,j)=a(i,j)-lnq1*a(i+1,j)-lnq2*a(i+2,j)+knq1*a(i-1,j)
enddo
enddo
end subroutine mf21b2i
subroutine mf22b2i(a,lnq1,lnq2,knq1,knq2, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2,knq1,knq2
integer :: i,j
do j=jts,jte
do i=ite,its,-1
a(i,j)=a(i,j)-lnq1*a(i+1,j)-lnq2*a(i+2,j) &
+knq1*a(i-1,j)+knq2*a(i-2,j)
enddo
enddo
end subroutine mf22b2i
subroutine ma11b2i(a,p,lnq1,knq1,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(1, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,knq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In mal1b2i; fnq exceeds tile width'
do j=jts,jte
p(1,j)=0.
enddo
do j=jts,jte
do i=its-1+fnq,its,-1
p(1,j)=a(i,j)-lnq1*p(1,j)+knq1*a(i-1,j)
enddo
enddo
end subroutine ma11b2i
subroutine ma12b2i(a,p,lnq1,knq1,knq2,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(1, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,knq1,knq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In mal2b2i; fnq exceeds tile width'
do j=jts,jte
p(1,j)=0.
enddo
do j=jts,jte
do i=its-1+fnq,its,-1
p(1,j)=a(i,j)-lnq1*p(1, j) &
+knq1*a(i-1,j)+knq2*a(i-2,j)
enddo
enddo
end subroutine ma12b2i
subroutine ma21b2i(a,p,lnq1,lnq2,knq1,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(2, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2,knq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In ma21b2i; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ma21b2i; fnq is not even'
do j=jts,jte
p(1,j)=0.
p(2,j)=0.
enddo
do j=jts,jte
do i=its-1+fnq,its,-2
p(2,j)=a(i, j)-lnq1*p(1, j)-lnq2*p(2,j) &
+knq1*a(i-1,j)
p(1,j)=a(i-1,j)-lnq1*p(2, j)-lnq2*p(1,j) &
+knq1*a(i-2,j)
enddo
enddo
end subroutine ma21b2i
subroutine ma22b2i(a,p,lnq1,lnq2,knq1,knq2,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(2, jms:jme), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2,knq1,knq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > ite+1-its)stop 'In ma22b2i; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ma22b2i; fnq is not even'
do j=jts,jte
p(1,j)=0.
p(2,j)=0.
enddo
do j=jts,jte
do i=its-1+fnq,its,-2
p(2,j)=a(i, j)-lnq1*p(1, j)-lnq2*p(2, j) &
+knq1*a(i-1,j)+knq2*a(i-2,j)
p(1,j)=a(i-1,j)-lnq1*p(2, j)-lnq2*p(1, j) &
+knq1*a(i-2,j)+knq2*a(i-3,j)
enddo
enddo
end subroutine ma22b2i
!############# 2-dimensions, index-j active:
subroutine ud1c2j(c,d,bnd1, & 1
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: c
real, dimension(ims:ime , jms:jme), &
INTENT(OUT ) :: d
real, INTENT(IN ) :: bnd1
integer :: i,j
do j=jts,jte
do i=its,ite
d(i,j)=bnd1*(c(i,j+1)-c(i,j-1))
enddo
enddo
end subroutine ud1c2j
subroutine ud2c2j(c,d,bnd1,bnd2, & 1
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: c
real, dimension(ims:ime , jms:jme), &
INTENT(OUT ) :: d
real, INTENT(IN ) :: bnd1,bnd2
integer :: i,j
do j=jts,jte
do i=its,ite
d(i,j)=bnd1*(c(i,j+1)-c(i,j-1))+bnd2*(c(i,j+2)-c(i,j-2))
enddo
enddo
end subroutine ud2c2j
subroutine sa1a2j(a,bnm1, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: bnm1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=bnm1*(a(i,j)+a(i,j+1))
enddo
enddo
end subroutine sa1a2j
subroutine sq1a2j(a,bnq, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: bnq
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j-1)+bnq*a(i,j)
enddo
enddo
end subroutine sq1a2j
subroutine sd1b2j(a,bnqi, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: bnqi
integer :: i,j
do j=jte,jts-1,-1
do i=its,ite
a(i,j)=bnqi*(a(i,j)-a(i,j-1))
enddo
enddo
end subroutine sd1b2j
subroutine ef1a2j(a,knq1, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: knq1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)+knq1*a(i,j+1)
enddo
enddo
end subroutine ef1a2j
subroutine ef2a2j(a,knq1,knq2, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: knq1,knq2
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)+knq1*a(i,j+1)+knq2*a(i,j+2)
enddo
enddo
end subroutine ef2a2j
subroutine ef1b2j(a,knq1, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: knq1
integer :: i,j
do j=jte,jts-1,-1
do i=its,ite
a(i,j)=a(i,j)+knq1*a(i,j-1)
enddo
enddo
end subroutine ef1b2j
subroutine ef2b2j(a,knq1,knq2, & 4
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: knq1,knq2
integer :: i,j
do j=jte,jts,-1
do i=its,ite
a(i,j)=a(i,j)+knq1*a(i,j-1)+knq2*a(i,j-2)
enddo
enddo
end subroutine ef2b2j
subroutine rf1a2j(a,lnq1, & 8
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j-1)
enddo
enddo
end subroutine rf1a2j
subroutine rf2a2j(a,lnq1,lnq2, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j-1)-lnq2*a(i,j-2)
enddo
enddo
end subroutine rf2a2j
subroutine rf1b2j(a,lnq1, & 8
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1
integer :: i,j
do j=jte,jts,-1
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j+1)
enddo
enddo
end subroutine rf1b2j
subroutine rf2b2j(a,lnq1,lnq2, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2
integer :: i,j
do j=jte,jts,-1
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j+1)-lnq2*a(i,j+2)
enddo
enddo
end subroutine rf2b2j
subroutine ra1a2j(a,p,lnq1,fnq, & 8
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 1), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In rala2j; fnq exceeds tile width'
do i=its,ite
p(i,1)=0.
enddo
do j=jte+1-fnq,jte,2
do i=its,ite
p(i,1)=a(i,j)-lnq1*p(i,1)
enddo
enddo
end subroutine ra1a2j
subroutine ra2a2j(a,p,lnq1,lnq2,fnq, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 2), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In ra2a2j; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ra2a2j; fnq is not even'
do j=1,2
do i=its,ite
p(i,j)=0.
enddo
enddo
do j=jte+1-fnq,jte,2
do i=its,ite
p(i,2)=a(i,j )-lnq1*p(i,1)-lnq2*p(i,2)
p(i,1)=a(i,j+1)-lnq1*p(i,2)-lnq2*p(i,1)
enddo
enddo
end subroutine ra2a2j
subroutine ra1b2j(a,p,lnq1,fnq, & 8
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 1), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In ralb2j; fnq exceeds tile width'
do i=its,ite
p(i,1)=0.
enddo
do j=jts-1+fnq,jts,-1
do i=its,ite
p(i,1)=a(i,j)-lnq1*p(i,1)
enddo
enddo
end subroutine ra1b2j
subroutine ra2b2j(a,p,lnq1,lnq2,fnq, & 6
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 2), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In ra2b2j; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ra2b2j; fnq is not even'
do j=1,2
do i=its,ite
p(i,j)=0.
enddo
enddo
do j=jts-1+fnq,jts,-2
do i=its,ite
p(i,2)=a(i,j )-lnq1*p(i,1)-lnq2*p(i,2)
p(i,1)=a(i,j-1)-lnq1*p(i,2)-lnq2*p(i,1)
enddo
enddo
end subroutine ra2b2j
subroutine mf11a2j(a,lnq1,knq1, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,knq1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j-1)+knq1*a(i,j+1)
enddo
enddo
end subroutine mf11a2j
subroutine mf12a2j(a,lnq1,knq1,knq2, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,knq1,knq2
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j-1) &
+knq1*a(i,j+1)+knq2*a(i,j+2)
enddo
enddo
end subroutine mf12a2j
subroutine mf21a2j(a,lnq1,lnq2,knq1, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2,knq1
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j-1)-lnq2*a(i,j-2)+knq1*a(i,j+1)
enddo
enddo
end subroutine mf21a2j
subroutine mf22a2j(a,lnq1,lnq2,knq1,knq2, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2,knq1,knq2
integer :: i,j
do j=jts,jte
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j-1)-lnq2*a(i,j-2) &
+knq1*a(i,j+1)+knq2*a(i,j+2)
enddo
enddo
end subroutine mf22a2j
subroutine ma11a2j(a,p,lnq1,knq1,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 1), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,knq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In mal1a2j; fnq exceeds tile width'
do i=its,ite
p(i,1)=0.
enddo
do j=jte+1-fnq,jte
do i=its,ite
p(i,1)=a(i,j)-lnq1*p(i,1)+knq1*a(i,j+1)
enddo
enddo
end subroutine ma11a2j
subroutine ma12a2j(a,p,lnq1,knq1,knq2,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 1), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,knq1,knq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In ma12a2j; fnq exceeds tile width'
do i=its,ite
p(i,1)=0.
enddo
do j=jte+1-fnq,jte
do i=its,ite
p(i,1)=a(i,j)-lnq1*p(i,1) &
+knq1*a(i,j+1)+knq2*a(i,j+2)
enddo
enddo
end subroutine ma12a2j
subroutine ma21a2j(a,p,lnq1,lnq2,knq1,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 2), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2,knq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In ma21a2j; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ma21a2j; fnq is not even'
do j=1,2
do i=its,ite
p(i,j)=0.
enddo
enddo
do j=jte+1-fnq,jte,2
do i=its,ite
p(i,2)=a(i,j )-lnq1*p(i,1)-lnq2*p(i,2) &
+knq1*a(i,j+1)
p(i,1)=a(i,j+1)-lnq1*p(i,2)-lnq2*p(i,1) &
+knq1*a(i,j+2)
enddo
enddo
end subroutine ma21a2j
subroutine ma22a2j(a,p,lnq1,lnq2,knq1,knq2,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 2), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2,knq1,knq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In ma22a2j; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ma22a2j; fnq is not even'
do j=1,2
do i=its,ite
p(i,j)=0.
enddo
enddo
do j=jte+1-fnq,jte,2
do i=its,ite
p(i,2)=a(i,j )-lnq1*p(i,1 )-lnq2*p(i,2 ) &
+knq1*a(i,j+1)+knq2*a(i,j+2)
p(i,1)=a(i,j+1)-lnq1*p(i,2 )-lnq2*p(i,1 ) &
+knq1*a(i,j+2)+knq2*a(i,j+3)
enddo
enddo
end subroutine ma22a2j
subroutine mf11b2j(a,lnq1,knq1, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,knq1
integer :: i,j
do j=jte,jts,-1
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j+1)+knq1*a(i,j-1)
enddo
enddo
end subroutine mf11b2j
subroutine mf12b2j(a,lnq1,knq1,knq2, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,knq1,knq2
integer :: i,j
do j=jte,jts,-1
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j+1) &
+knq1*a(i,j-1)+knq2*a(i,j-2)
enddo
enddo
end subroutine mf12b2j
subroutine mf21b2j(a,lnq1,lnq2,knq1, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2,knq1
integer :: i,j
do j=jte,jts,-1
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j+1)-lnq2*a(i,j+2)+knq1*a(i,j-1)
enddo
enddo
end subroutine mf21b2j
subroutine mf22b2j(a,lnq1,lnq2,knq1,knq2, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(INOUT) :: a
real, INTENT(IN ) :: lnq1,lnq2,knq1,knq2
integer :: i,j
do j=jte,jts,-1
do i=its,ite
a(i,j)=a(i,j)-lnq1*a(i,j+1)-lnq2*a(i,j+2) &
+knq1*a(i,j-1)+knq2*a(i,j-2)
enddo
enddo
end subroutine mf22b2j
subroutine ma11b2j(a,p,lnq1,knq1,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 1),&
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,knq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In mal1b2j; fnq exceeds tile width'
do i=its,ite
p(i,1)=0.
enddo
do j=jts-1+fnq,jts,-1
do i=its,ite
p(i,1)=a(i,j)-lnq1*p(i,1)+knq1*a(i,j-1)
enddo
enddo
end subroutine ma11b2j
subroutine ma12b2j(a,p,lnq1,knq1,knq2,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 1), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,knq1,knq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In mal2b2j; fnq exceeds tile width'
do i=its,ite
p(i,1)=0.
enddo
do j=jts-1+fnq,jts,-1
do i=its,ite
p(i,1)=a(i,j)-lnq1*p(i,1 ) &
+knq1*a(i,j-1)+knq2*a(i,j-2)
enddo
enddo
end subroutine ma12b2j
subroutine ma21b2j(a,p,lnq1,lnq2,knq1,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime , 2), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2,knq1
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In ma21b2j; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ma21b2j; fnq is not even'
do j=1,2
do i=its,ite
p(i,j)=0.
enddo
enddo
do j=jts-1+fnq,jts,-2
do i=its,ite
p(i,2)=a(i,j )-lnq1*p(i,1 )-lnq2*p(i,2) &
+knq1*a(i,j-1)
p(i,1)=a(i,j-1)-lnq1*p(i,2 )-lnq2*p(i,1) &
+knq1*a(i,j-2)
enddo
enddo
end subroutine ma21b2j
subroutine ma22b2j(a,p,lnq1,lnq2,knq1,knq2,fnq, &
ids,ide, jds,jde, &
ims,ime, jms,jme, &
its,ite, jts,jte)
implicit none
INTEGER, INTENT(IN ) :: ids,ide, jds,jde
INTEGER, INTENT(IN ) :: ims,ime, jms,jme
INTEGER, INTENT(IN ) :: its,ite, jts,jte
real, dimension(ims:ime , jms:jme), &
INTENT(IN ) :: a
real, dimension(ims:ime, 2), &
INTENT(INOUT) :: p
real, INTENT(IN ) :: lnq1,lnq2,knq1,knq2
integer, INTENT(IN ) :: fnq
integer :: i,j
if(fnq > jte+1-jts)stop 'In ma22b2j; fnq exceeds tile width'
if(mod(fnq,2)/=0) stop 'In ma22b2j; fnq is not even'
do j=1,2
do i=its,ite
p(i,j)=0.
enddo
enddo
do j=jts-1+fnq,jts,-2
do i=its,ite
p(i,2)=a(i,j )-lnq1*p(i,1 )-lnq2*p(i,2 ) &
+knq1*a(i,j-1)+knq2*a(i,j-2)
p(i,1)=a(i,j-1)-lnq1*p(i,2 )-lnq2*p(i,1 ) &
+knq1*a(i,j-2)+knq2*a(i,j-3)
enddo
enddo
end subroutine ma22b2j