trilinear_interp_var Function

public function trilinear_interp_var(var, cell, xgrid, ygrid, zgrid, x, y, z)

Uses

  • proc~~trilinear_interp_var~~UsesGraph proc~trilinear_interp_var trilinear_interp_var decomp_2d decomp_2d proc~trilinear_interp_var->decomp_2d module~modglobal modglobal proc~trilinear_interp_var->module~modglobal

Arguments

Type IntentOptional Attributes Name
real, intent(in) :: var(ib-ih:ie+ih,jb-jh:je+jh,kb-kh:kb+kh)
integer, intent(in) :: cell(3)
real, intent(in), dimension(ib:itot+ih) :: xgrid
real, intent(in), dimension(jb:jtot+jh) :: ygrid
real, intent(in), dimension(kb:ktot+kh) :: zgrid
real, intent(in) :: x
real, intent(in) :: y
real, intent(in) :: z

Return Value real


Calls

proc~~trilinear_interp_var~~CallsGraph proc~trilinear_interp_var trilinear_interp_var proc~eval_corners eval_corners proc~trilinear_interp_var->proc~eval_corners proc~trilinear_interp trilinear_interp proc~trilinear_interp_var->proc~trilinear_interp zstart zstart proc~trilinear_interp_var->zstart

Called by

proc~~trilinear_interp_var~~CalledByGraph proc~trilinear_interp_var trilinear_interp_var proc~wallfunheat wallfunheat proc~wallfunheat->proc~trilinear_interp_var proc~wallfunmom wallfunmom proc~wallfunmom->proc~trilinear_interp_var proc~ibmwallfun ibmwallfun proc~ibmwallfun->proc~wallfunheat proc~ibmwallfun->proc~wallfunmom program~dalesurban DALESURBAN program~dalesurban->proc~ibmwallfun

Source Code

   real function trilinear_interp_var(var, cell, xgrid, ygrid, zgrid, x, y, z)
     use modglobal, only : ib, ie, ih, jb, je, jh, kb, ke, kh, itot, jtot, ktot
     use decomp_2d, only : zstart
     implicit none
     real, intent(in)    :: var(ib-ih:ie+ih,jb-jh:je+jh,kb-kh:kb+kh)
     integer, intent(in) :: cell(3) ! GLOBAL indices of cell containing the point
     real, intent(in), dimension(ib:itot+ih) :: xgrid
     real, intent(in), dimension(jb:jtot+jh) :: ygrid
     real, intent(in), dimension(kb:ktot+kh) :: zgrid
     real,    intent(in) :: x, y, z ! location of point to interpolate at
     real, dimension(8)  :: corners(8)
     real :: x0, y0, z0, x1, y1, z1
     integer :: i, j, k

     i = cell(1) - zstart(1) + 1
     j = cell(2) - zstart(2) + 1
     k = cell(3) - zstart(3) + 1
     if ((i < ib-1) .or. (i > ie+1) .or. (j < jb-1) .or. (j > je+1)) then
       write(*,*) "problem in trilinear_interp_var", i, j, k
       stop 1
     end if
     corners = eval_corners(var, i, j, k)

     x0 = xgrid(cell(1))
     y0 = ygrid(cell(2))
     z0 = zgrid(cell(3))
     x1 = xgrid(cell(1)+1)
     y1 = ygrid(cell(2)+1)
     z1 = zgrid(cell(3)+1)

     trilinear_interp_var = trilinear_interp(x, y, z, x0, y0, z0, x1, y1, z1, corners)

   end function trilinear_interp_var