cubt4u Subroutine

private subroutine cubt4u(nx, ny, nz, nt, p, mx, my, mz, mt, mxmy, mxmymz, q, intpol, lt, dtm, dt, dtp, dtpp, ptm, pt, ptp, ptpp, 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, lnmt, isubx, jsuby, ksubz, lsubt)

cubically interpolate in t

Arguments

Type IntentOptional Attributes Name
integer :: nx
integer :: ny
integer :: nz
integer :: nt
real(kind=wp) :: p(nx,ny,nz,nt)
integer :: mx
integer :: my
integer :: mz
integer :: mt
integer :: mxmy
integer :: mxmymz
real(kind=wp) :: q(mxmymz,mt)
integer :: intpol(4)
integer :: lt(mt)
real(kind=wp) :: dtm(mt)
real(kind=wp) :: dt(mt)
real(kind=wp) :: dtp(mt)
real(kind=wp) :: dtpp(mt)
real(kind=wp) :: ptm(mxmymz)
real(kind=wp) :: pt(mxmymz)
real(kind=wp) :: ptp(mxmymz)
real(kind=wp) :: ptpp(mxmymz)
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 :: lnmt
integer :: isubx
integer :: jsuby
integer :: ksubz
integer :: lsubt

Calls

proc~~cubt4u~~CallsGraph proc~cubt4u regridpack_module::cubt4u proc~cubt3u regridpack_module::cubt3u proc~cubt4u->proc~cubt3u proc~lint3u regridpack_module::lint3u proc~cubt4u->proc~lint3u proc~cubt2u regridpack_module::cubt2u proc~cubt3u->proc~cubt2u proc~lint2u regridpack_module::lint2u proc~cubt3u->proc~lint2u proc~lint3u->proc~cubt2u proc~lint3u->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~~cubt4u~~CalledByGraph proc~cubt4u regridpack_module::cubt4u proc~rgrd4u regridpack_module::rgrd4u proc~rgrd4u->proc~cubt4u interface~regrid regridpack_module::regrid interface~regrid->proc~rgrd4u

Source Code

    subroutine cubt4u(nx,ny,nz,nt,p,mx,my,mz,mt,mxmy,mxmymz,q,intpol,&
                        lt,dtm,dt,dtp,dtpp,ptm,pt,ptp,ptpp,&
                        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,lnmt,isubx,jsuby,ksubz,lsubt)

    implicit none

    integer  :: nx,ny,nz,nt,mx,my,mz,mt
    integer  :: inmx,jnmy,knmz,lnmt,isubx,jsuby,ksubz,lsubt
    integer  :: mxmy,mxmymz
    real(wp) :: p(nx,ny,nz,nt),q(mxmymz,mt)
    real(wp) :: ptm(mxmymz),pt(mxmymz),ptp(mxmymz),ptpp(mxmymz)
    real(wp) :: dtm(mt),dt(mt),dtp(mt),dtpp(mt)
    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)
    integer  :: lt(mt),kz(mz),jy(my),ix(mx),intpol(4)
    integer  :: l,ll,iijjkk,lsave

    if (intpol(3) == 1) then

        ! linear in z

        if (lsubt == 1) then
            ! mt grid is subset of nt grid
            do ll=1,mt
                l = lnmt*(ll-1)+1
                call lint3u(nx,ny,nz,p(1,1,1,l),mx,my,mxmy,mz,q(1,ll),intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
            end do
            return
        end if
        lsave = -3
        do ll=1,mt
            l = lt(ll)
            if (l==lsave) then
                ! l pointer has not moved since last pass (no updates or interpolation)
            else if (l==lsave+1) then
                ! update l-1,l,l+1 and interpolate l+2
                do iijjkk=1,mxmymz
                    ptm(iijjkk) = pt(iijjkk)
                    pt(iijjkk) = ptp(iijjkk)
                    ptp(iijjkk) = ptpp(iijjkk)
                end do
                call lint3u(nx,ny,nz,p(1,1,1,l+2),mx,my,mxmy,mz,ptpp,intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
            else if (l==lsave+2) then
                ! update l-1,l and interpolate l+1,l+2
                do iijjkk=1,mxmymz
                    ptm(iijjkk) = ptp(iijjkk)
                    pt(iijjkk) = ptpp(iijjkk)
                end do
                call lint3u(nx,ny,nz,p(1,1,1,l+1),mx,my,mxmy,mz,ptp,intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
                call lint3u(nx,ny,nz,p(1,1,1,l+2),mx,my,mxmy,mz,ptpp,intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
            else if (l==lsave+3) then
                ! update l-1 and interpolate l,l+1,l+2
                do iijjkk=1,mxmymz
                    ptm(iijjkk) = ptpp(iijjkk)
                end do
                call lint3u(nx,ny,nz,p(1,1,1,l),mx,my,mxmy,mz,pt,intpol,kz,dz,&
                            pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
                call lint3u(nx,ny,nz,p(1,1,1,l+1),mx,my,mxmy,mz,ptp,intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
                call lint3u(nx,ny,nz,p(1,1,1,l+2),mx,my,mxmy,mz,ptpp,intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
            else
                ! interpolate all four l-1,l,l+1,l+2
                call lint3u(nx,ny,nz,p(1,1,1,l-1),mx,my,mxmy,mz,ptm,intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
                call lint3u(nx,ny,nz,p(1,1,1,l),mx,my,mxmy,mz,pt,intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
                call lint3u(nx,ny,nz,p(1,1,1,l+1),mx,my,mxmy,mz,ptp,intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
                call lint3u(nx,ny,nz,p(1,1,1,l+2),mx,my,mxmy,mz,ptpp,intpol,kz,&
                            dz,pk,pkp,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp,&
                            inmx,jnmy,knmz,isubx,jsuby,ksubz)
            end if

            ! save l pointer for next pass

            lsave = l

            ! cubically interpolate q(ii,jj,kk,ll) from ptm,pt,ptp,ptpp in t direction

            do iijjkk=1,mxmymz
                q(iijjkk,ll) = dtm(ll)*ptm(iijjkk) + dt(ll)*pt(iijjkk) + dtp(ll)*ptp(iijjkk) + dtpp(ll)*ptpp(iijjkk)
            end do
        end do

    else

        ! cubic in z

        if (lsubt == 1) then

            ! mt grid is subset of nt grid

            do ll=1,mt
                l = lnmt*(ll-1)+1
                call cubt3u(nx,ny,nz,p(1,1,1,l),mx,my,mxmy,mz,q(1,ll),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)
            end do
            return
        end if
        lsave = -3
        do ll=1,mt
            l = lt(ll)
            if (l==lsave) then

                ! l pointer has not moved since last pass (no updates or interpolation)

            else if (l==lsave+1) then

                ! update l-1,l,l+1 and interpolate l+2

                do iijjkk=1,mxmymz
                    ptm(iijjkk) = pt(iijjkk)
                    pt(iijjkk) = ptp(iijjkk)
                    ptp(iijjkk) = ptpp(iijjkk)
                end do
                call cubt3u(nx,ny,nz,p(1,1,1,l+2),mx,my,mxmy,mz,ptpp,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)
            else if (l==lsave+2) then

                ! update l-1,l and interpolate l+1,l+2

                do iijjkk=1,mxmymz
                    ptm(iijjkk) = ptp(iijjkk)
                    pt(iijjkk) = ptpp(iijjkk)
                end do
                call cubt3u(nx,ny,nz,p(1,1,1,l+1),mx,my,mxmy,mz,ptp,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)
                call cubt3u(nx,ny,nz,p(1,1,1,l+2),mx,my,mxmy,mz,ptpp,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)
            else if (l==lsave+3) then

                ! update l-1 and interpolate l,l+1,l+2

                do iijjkk=1,mxmymz
                    ptm(iijjkk) = ptpp(iijjkk)
                end do
                call cubt3u(nx,ny,nz,p(1,1,1,l),mx,my,mxmy,mz,pt,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)
                call cubt3u(nx,ny,nz,p(1,1,1,l+1),mx,my,mxmy,mz,ptp,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)
                call cubt3u(nx,ny,nz,p(1,1,1,l+2),mx,my,mxmy,mz,ptpp,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)
            else

                ! interpolate all four l-1,l,l+1,l+2

                call cubt3u(nx,ny,nz,p(1,1,1,l-1),mx,my,mxmy,mz,ptm,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)
                call cubt3u(nx,ny,nz,p(1,1,1,l),mx,my,mxmy,mz,pt,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)
                call cubt3u(nx,ny,nz,p(1,1,1,l+1),mx,my,mxmy,mz,ptp,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)
                call cubt3u(nx,ny,nz,p(1,1,1,l+2),mx,my,mxmy,mz,ptpp,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)
            end if

            ! save l pointer for next pass
            lsave = l

            ! cubically interpolate q(ii,jj,kk,ll) from ptm,pt,ptp,ptpp in t direction
            do iijjkk=1,mxmymz
                q(iijjkk,ll) = dtm(ll)*ptm(iijjkk) + dt(ll)*pt(iijjkk) + dtp(ll)*ptp(iijjkk) + dtpp(ll)*ptpp(iijjkk)
            end do
        end do

    end if

    end subroutine cubt4u