cubt2u Subroutine

private subroutine cubt2u(nx, ny, p, mx, my, q, intpol, jy, dym, dy, dyp, dypp, pjm, pj, pjp, pjpp, ix, dxm, dx, dxp, dxpp, inmx, jnmy, isubx, jsuby)

cubically interpolate p onto q in y

Arguments

Type IntentOptional Attributes Name
integer :: nx
integer :: ny
real(kind=wp) :: p(nx,ny)
integer :: mx
integer :: my
real(kind=wp) :: q(mx,my)
integer :: intpol(2)
integer :: jy(my)
real(kind=wp) :: dym(my)
real(kind=wp) :: dy(my)
real(kind=wp) :: dyp(my)
real(kind=wp) :: dypp(my)
real(kind=wp) :: pjm(mx)
real(kind=wp) :: pj(mx)
real(kind=wp) :: pjp(mx)
real(kind=wp) :: pjpp(mx)
integer :: ix(mx)
real(kind=wp) :: dxm(mx)
real(kind=wp) :: dx(mx)
real(kind=wp) :: dxp(mx)
real(kind=wp) :: dxpp(mx)
integer :: inmx
integer :: jnmy
integer :: isubx
integer :: jsuby

Calls

proc~~cubt2u~~CallsGraph proc~cubt2u regridpack_module::cubt2u proc~cubt1u regridpack_module::cubt1u proc~cubt2u->proc~cubt1u proc~lint1u regridpack_module::lint1u proc~cubt2u->proc~lint1u

Called by

proc~~cubt2u~~CalledByGraph proc~cubt2u regridpack_module::cubt2u proc~cubt3u regridpack_module::cubt3u proc~cubt3u->proc~cubt2u proc~lint3u regridpack_module::lint3u proc~lint3u->proc~cubt2u proc~rgrd2u regridpack_module::rgrd2u proc~rgrd2u->proc~cubt2u interface~regrid regridpack_module::regrid interface~regrid->proc~rgrd2u proc~rgrd3u regridpack_module::rgrd3u interface~regrid->proc~rgrd3u proc~rgrd4u regridpack_module::rgrd4u interface~regrid->proc~rgrd4u proc~cubt4u regridpack_module::cubt4u proc~cubt4u->proc~cubt3u proc~cubt4u->proc~lint3u proc~lint4u regridpack_module::lint4u proc~lint4u->proc~cubt3u proc~lint4u->proc~lint3u proc~rgrd3u->proc~cubt3u proc~rgrd3u->proc~lint3u proc~rgrd4u->proc~cubt4u proc~rgrd4u->proc~lint4u

Source Code

    subroutine cubt2u(nx,ny,p,mx,my,q,intpol,jy,dym,dy,dyp,dypp,&
                      pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,&
                      isubx,jsuby)

    implicit none

    integer  :: nx,ny,mx,my,intpol(2),jy(my),ix(mx),inmx,jnmy,isubx,jsuby
    integer  :: j,jj,ii,jsave
    real(wp) :: p(nx,ny),q(mx,my)
    real(wp) :: dym(my),dy(my),dyp(my),dypp(my)
    real(wp) :: pjm(mx),pj(mx),pjp(mx),pjpp(mx)
    real(wp) :: dxm(mx),dx(mx),dxp(mx),dxpp(mx)

    if (intpol(1) == 1) then

        ! linear in x

        if (jsuby == 1) then
            ! my grid is subset of ny grid
            do jj=1,my
                j = jnmy*(jj-1)+1
                call lint1u(nx,p(1,j),mx,q(1,jj),ix,dx,inmx,isubx)
            end do
            return
        end if

        jsave = -3
        do jj=1,my
            j = jy(jj)

            ! load pjm,pj,pjp,pjpp

            if (j==jsave) then
                ! no updates or x interpolation necessary
            else if (j==jsave+1) then
                do ii=1,mx
                    pjm(ii) = pj(ii)
                    pj(ii) = pjp(ii)
                    pjp(ii) = pjpp(ii)
                end do
                call lint1u(nx,p(1,j+2),mx,pjpp,ix,dx,inmx,isubx)
            else if (j==jsave+2) then
                do ii=1,mx
                    pjm(ii) = pjp(ii)
                    pj(ii) = pjpp(ii)
                end do
                call lint1u(nx,p(1,j+1),mx,pjp,ix,dx,inmx,isubx)
                call lint1u(nx,p(1,j+2),mx,pjpp,ix,dx,inmx,isubx)
            else if (j==jsave+3) then
                do ii=1,mx
                    pjm(ii) = pjpp(ii)
                end do
                call lint1u(nx,p(1,j),mx,pj,ix,dx,inmx,isubx)
                call lint1u(nx,p(1,j+1),mx,pjp,ix,dx,inmx,isubx)
                call lint1u(nx,p(1,j+2),mx,pjpp,ix,dx,inmx,isubx)
            else
                ! load all four (no updates)
                call lint1u(nx,p(1,j-1),mx,pjm,ix,dx,inmx,isubx)
                call lint1u(nx,p(1,j),mx,pj,ix,dx,inmx,isubx)
                call lint1u(nx,p(1,j+1),mx,pjp,ix,dx,inmx,isubx)
                call lint1u(nx,p(1,j+2),mx,pjpp,ix,dx,inmx,isubx)
            end if
            ! update pointer
            jsave = j
            do ii=1,mx
                q(ii,jj) = dym(jj)*pjm(ii) + dy(jj)*pj(ii) + dyp(jj)*pjp(ii) + dypp(jj)*pjpp(ii)
            end do
        end do

    else

        ! cubic in x

        if (jsuby == 1) then
            ! my grid is subset of ny grid
            do jj=1,my
                j = jnmy*(jj-1)+1
                call cubt1u(nx,p(1,j),mx,q(1,jj),ix,dxm,dx,dxp,dxpp,inmx,isubx)
            end do
            return
        end if

        jsave = -3
        do jj=1,my
            j = jy(jj)

            ! load pjm,pj,pjp,pjpp

            if (j==jsave) then
                ! no updates or x interpolation necessary
            else if (j==jsave+1) then
                do ii=1,mx
                    pjm(ii) = pj(ii)
                    pj(ii) = pjp(ii)
                    pjp(ii) = pjpp(ii)
                end do
                call cubt1u(nx,p(1,j+2),mx,pjpp,ix,dxm,dx,dxp,dxpp,inmx,isubx)
            else if (j==jsave+2) then
                do ii=1,mx
                    pjm(ii) = pjp(ii)
                    pj(ii) = pjpp(ii)
                end do
                call cubt1u(nx,p(1,j+1),mx,pjp,ix,dxm,dx,dxp,dxpp,inmx,isubx)
                call cubt1u(nx,p(1,j+2),mx,pjpp,ix,dxm,dx,dxp,dxpp,inmx,isubx)
            else if (j==jsave+3) then
                do ii=1,mx
                    pjm(ii) = pjpp(ii)
                end do
                call cubt1u(nx,p(1,j),mx,pj,ix,dxm,dx,dxp,dxpp,inmx,isubx)
                call cubt1u(nx,p(1,j+1),mx,pjp,ix,dxm,dx,dxp,dxpp,inmx,isubx)
                call cubt1u(nx,p(1,j+2),mx,pjpp,ix,dxm,dx,dxp,dxpp,inmx,isubx)
            else
                ! load all four (no updates)
                call cubt1u(nx,p(1,j-1),mx,pjm,ix,dxm,dx,dxp,dxpp,inmx,isubx)
                call cubt1u(nx,p(1,j),mx,pj,ix,dxm,dx,dxp,dxpp,inmx,isubx)
                call cubt1u(nx,p(1,j+1),mx,pjp,ix,dxm,dx,dxp,dxpp,inmx,isubx)
                call cubt1u(nx,p(1,j+2),mx,pjpp,ix,dxm,dx,dxp,dxpp,inmx,isubx)
            end if
            ! update pointer
            jsave = j
            do ii=1,mx
                q(ii,jj) = dym(jj)*pjm(ii) + dy(jj)*pj(ii) + dyp(jj)*pjp(ii) + dypp(jj)*pjpp(ii)
            end do
        end do
        return

    end if

    end subroutine cubt2u