cubt1u Subroutine

private subroutine cubt1u(nx, p, mx, q, ix, dxm, dx, dxp, dxpp, inmx, isubx)

Arguments

Type IntentOptional Attributes Name
integer :: nx
real(kind=wp) :: p(nx)
integer :: mx
real(kind=wp) :: q(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 :: isubx

Called by

proc~~cubt1u~~CalledByGraph proc~cubt1u regridpack_module::cubt1u proc~cubt2u regridpack_module::cubt2u proc~cubt2u->proc~cubt1u proc~lint2u regridpack_module::lint2u proc~lint2u->proc~cubt1u proc~rgrd1u regridpack_module::rgrd1u proc~rgrd1u->proc~cubt1u interface~regrid regridpack_module::regrid interface~regrid->proc~rgrd1u proc~rgrd2u regridpack_module::rgrd2u interface~regrid->proc~rgrd2u proc~rgrd3u regridpack_module::rgrd3u interface~regrid->proc~rgrd3u proc~rgrd4u regridpack_module::rgrd4u interface~regrid->proc~rgrd4u proc~cubt3u regridpack_module::cubt3u proc~cubt3u->proc~cubt2u proc~cubt3u->proc~lint2u proc~lint3u regridpack_module::lint3u proc~lint3u->proc~cubt2u proc~lint3u->proc~lint2u proc~rgrd2u->proc~cubt2u proc~rgrd2u->proc~lint2u 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 cubt1u(nx,p,mx,q,ix,dxm,dx,dxp,dxpp,inmx,isubx)

    implicit none

    integer  :: nx,mx,ix(mx),inmx,isubx,i,ii
    real(wp) :: p(nx),q(mx)
    real(wp) :: dxm(mx),dx(mx),dxp(mx),dxpp(mx)

    if (isubx == 1) then
        ! mx grid is subset of nx grid so q can be set directly
        do ii=1,mx
            i = inmx*(ii-1)+1
            q(ii) = p(i)
        end do
    else
        ! cubically interpolate on uniform grid
        do ii=1,mx
            i = ix(ii)
            q(ii)=(dxm(ii)*p(i-1)+dx(ii)*p(i)+dxp(ii)*p(i+1)+dxpp(ii)*p(i+2))
        end do
    end if

    end subroutine cubt1u