cubt3u Subroutine

private subroutine cubt3u(nx, ny, nz, p, mx, my, mxmy, mz, q, intpol, kz, dzm, dz, dzp, dzpp, pkm, pk, pkp, pkpp, jy, dym, dy, dyp, dypp, pjm, pj, pjp, pjpp, ix, dxm, dx, dxp, dxpp, inmx, jnmy, knmz, isubx, jsuby, ksubz)

cubically interpolate in z

Arguments

Type IntentOptional Attributes Name
integer :: nx
integer :: ny
integer :: nz
real(kind=wp) :: p(nx,ny,nz)
integer :: mx
integer :: my
integer :: mxmy
integer :: mz
real(kind=wp) :: q(mxmy,mz)
integer :: intpol(3)
integer :: kz(mz)
real(kind=wp) :: dzm(mz)
real(kind=wp) :: dz(mz)
real(kind=wp) :: dzp(mz)
real(kind=wp) :: dzpp(mz)
real(kind=wp) :: pkm(mxmy)
real(kind=wp) :: pk(mxmy)
real(kind=wp) :: pkp(mxmy)
real(kind=wp) :: pkpp(mxmy)
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 :: knmz
integer :: isubx
integer :: jsuby
integer :: ksubz

Calls

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

Called by

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

Source Code

    subroutine cubt3u(nx,ny,nz,p,mx,my,mxmy,mz,q,intpol,&
                      kz,dzm,dz,dzp,dzpp,pkm,pk,pkp,pkpp,jy,dym,dy,dyp,dypp,pjm,pj,&
                      pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                      inmx,jnmy,knmz,isubx,jsuby,ksubz)

    implicit none

    integer  :: nx,ny,nz,mx,my,mz,mxmy,intpol(3),kz(mz),jy(my),ix(mx)
    integer  :: inmx,jnmy,knmz,isubx,jsuby,ksubz
    integer  :: kk,k,iijj,ksave
    real(wp) :: p(nx,ny,nz),q(mxmy,mz)
    real(wp) :: dzm(mz),dz(mz),dzp(mz),dzpp(mz)
    real(wp) :: pkm(mxmy),pk(mxmy),pkp(mxmy),pkpp(mxmy)
    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(2) == 1) then

        ! linear in y
        if (ksubz == 1) then
            ! mz grid is subset of nz grid
            do kk=1,mz
                k = knmz*(kk-1)+1
                call lint2u(nx,ny,p(1,1,k),mx,my,q(1,kk),intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
            end do
            return
        end if

        ! mz not a subgrid of nz

        ksave = -3
        do kk=1,mz
            k = kz(kk)
            if (k==ksave) then
                ! k pointer has not moved since last pass (no updates or interpolation)
            else if (k==ksave+1) then
                ! update k-1,k,k+1 and interpolate k+2
                do iijj=1,mxmy
                    pkm(iijj) = pk(iijj)
                    pk(iijj) = pkp(iijj)
                    pkp(iijj) = pkpp(iijj)
                end do
                call lint2u(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
            else if (k==ksave+2) then
                ! update k-1,k and interpolate k+1,k+2
                do iijj=1,mxmy
                    pkm(iijj) = pkp(iijj)
                    pk(iijj) = pkpp(iijj)
                end do
                call lint2u(nx,ny,p(1,1,k+1),mx,my,pkp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
                call lint2u(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
            else if (k==ksave+3) then
                ! update k-1 and interpolate k,k+1,k+2
                do iijj=1,mxmy
                    pkm(iijj) = pkpp(iijj)
                end do
                call lint2u(nx,ny,p(1,1,k),mx,my,pk,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
                call lint2u(nx,ny,p(1,1,k+1),mx,my,pkp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
                call lint2u(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
            else
                ! interpolate all four k-1,k,k+1,k+2
                call lint2u(nx,ny,p(1,1,k-1),mx,my,pkm,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
                call lint2u(nx,ny,p(1,1,k),mx,my,pk,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
                call lint2u(nx,ny,p(1,1,k+1),mx,my,pkp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
                call lint2u(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp,inmx,jnmy,isubx,jsuby)
            end if

            ! save k pointer for next pass

            ksave = k

            ! cubically interpolate q(ii,jj,kk) from pkm,pk,pkp,pkpp in z direction

            do iijj=1,mxmy
                q(iijj,kk) = dzm(kk)*pkm(iijj) + dz(kk)*pk(iijj) + dzp(kk)*pkp(iijj) + dzpp(kk)*pkpp(iijj)
            end do
        end do

    else

        ! cubic in y

        if (ksubz == 1) then
            ! mz grid is subset of nz grid
            do kk=1,mz
                k = knmz*(kk-1)+1
                call cubt2u(nx,ny,p(1,1,k),mx,my,q(1,kk),intpol,jy,dym,dy,&
                    dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                    inmx,jnmy,isubx,jsuby)

            end do
            return
        end if

        ksave = -3
        do kk=1,mz
            k = kz(kk)
            if (k==ksave) then
                ! k pointer has not moved since last pass (no updates or interpolation)
            else if (k==ksave+1) then
                ! update k-1,k,k+1 and interpolate k+2
                do iijj=1,mxmy
                    pkm(iijj) = pk(iijj)
                    pk(iijj) = pkp(iijj)
                    pkp(iijj) = pkpp(iijj)
                end do
                call cubt2u(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
            else if (k==ksave+2) then
                ! update k-1,k and interpolate k+1,k+2
                do iijj=1,mxmy
                    pkm(iijj) = pkp(iijj)
                    pk(iijj) = pkpp(iijj)
                end do
                call cubt2u(nx,ny,p(1,1,k+1),mx,my,pkp,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
                call cubt2u(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
            else if (k==ksave+3) then
                ! update k-1 and interpolate k,k+1,k+2
                do iijj=1,mxmy
                    pkm(iijj) = pkpp(iijj)
                end do
                call cubt2u(nx,ny,p(1,1,k),mx,my,pk,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
                call cubt2u(nx,ny,p(1,1,k+1),mx,my,pkp,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
                call cubt2u(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
            else
                ! interpolate all four k-1,k,k+1,k+2
                call cubt2u(nx,ny,p(1,1,k-1),mx,my,pkm,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
                call cubt2u(nx,ny,p(1,1,k),mx,my,pk,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
                call cubt2u(nx,ny,p(1,1,k+1),mx,my,pkp,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
                call cubt2u(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,&
                            jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,isubx,jsuby)
            end if

            ! save k pointer for next pass

            ksave = k

            ! cubically interpolate q(ii,jj,kk) from pkm,pk,pkp,pkpp in z direction

            do iijj=1,mxmy
                q(iijj,kk) = dzm(kk)*pkm(iijj) + dz(kk)*pk(iijj) + dzp(kk)*pkp(iijj) + dzpp(kk)*pkpp(iijj)
            end do
        end do

    end if

    end subroutine cubt3u