cubt3 Subroutine

private subroutine cubt3(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)

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)

Calls

proc~~cubt3~~CallsGraph proc~cubt3 regridpack_module::cubt3 proc~cubt2 regridpack_module::cubt2 proc~cubt3->proc~cubt2 proc~lint2 regridpack_module::lint2 proc~cubt3->proc~lint2 proc~cubt1 regridpack_module::cubt1 proc~cubt2->proc~cubt1 proc~lint1 regridpack_module::lint1 proc~cubt2->proc~lint1 proc~lint2->proc~cubt1 proc~lint2->proc~lint1

Called by

proc~~cubt3~~CalledByGraph proc~cubt3 regridpack_module::cubt3 proc~cubt4 regridpack_module::cubt4 proc~cubt4->proc~cubt3 proc~lint4 regridpack_module::lint4 proc~lint4->proc~cubt3 proc~rgrd3 regridpack_module::rgrd3 proc~rgrd3->proc~cubt3 interface~regrid regridpack_module::regrid interface~regrid->proc~rgrd3 proc~rgrd3_wrapper regridpack_module::rgrd3_wrapper interface~regrid->proc~rgrd3_wrapper proc~rgrd4 regridpack_module::rgrd4 interface~regrid->proc~rgrd4 proc~rgrd4_wrapper regridpack_module::rgrd4_wrapper interface~regrid->proc~rgrd4_wrapper proc~rgrd3_wrapper->proc~rgrd3 proc~rgrd4->proc~cubt4 proc~rgrd4->proc~lint4 proc~rgrd4_wrapper->proc~rgrd4

Source Code

    subroutine cubt3(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)

    implicit none

    integer  :: nx,ny,nz,mx,my,mxmy,mz,k,kk,ksave,iijj
    real(wp) :: p(nx,ny,nz),q(mxmy,mz)
    real(wp) :: pkm(mxmy),pk(mxmy),pkp(mxmy),pkpp(mxmy)
    real(wp) :: pjm(mx),pj(mx),pjp(mx),pjpp(mx)
    real(wp) :: dzm(mz),dz(mz),dzp(mz),dzpp(mz)
    real(wp) :: dym(my),dy(my),dyp(my),dypp(my)
    real(wp) :: dxm(mx),dx(mx),dxp(mx),dxpp(mx)
    integer  :: intpol(3),kz(mz),jy(my),ix(mx)

    if (intpol(2) == 1) then

        ! linear in y

        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 lint2(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
            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 lint2(nx,ny,p(1,1,k+1),mx,my,pkp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
                call lint2(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
            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 lint2(nx,ny,p(1,1,k),mx,my,pk,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
                call lint2(nx,ny,p(1,1,k+1),mx,my,pkp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
                call lint2(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
            else
                ! interpolate all four k-1,k,k+1,k+2
                call lint2(nx,ny,p(1,1,k-1),mx,my,pkm,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
                call lint2(nx,ny,p(1,1,k),mx,my,pk,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
                call lint2(nx,ny,p(1,1,k+1),mx,my,pkp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
                call lint2(nx,ny,p(1,1,k+2),mx,my,pkpp,intpol,jy,dy,pj,pjp,ix,dxm,dx,dxp,dxpp)
            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

        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 cubt2(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)
            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 cubt2(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)
                call cubt2(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)
            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 cubt2(nx,ny,p(1,1,k),mx,my,pk,intpol,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp)
                call cubt2(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)
                call cubt2(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)
            else
                ! interpolate all four k-1,k,k+1,k+2
                call cubt2(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)
                call cubt2(nx,ny,p(1,1,k),mx,my,pk,intpol,jy,dym,dy,dyp,dypp,pjm,pj,pjp,pjpp,ix,dxm,dx,dxp,dxpp)
                call cubt2(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)
                call cubt2(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)
            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 cubt3