Skip to content

Commit

Permalink
adding subgrid ratio to gradient calculation
Browse files Browse the repository at this point in the history
  • Loading branch information
jbeezley committed Jun 14, 2010
1 parent e3dd2bd commit 80d70fb
Showing 1 changed file with 28 additions and 18 deletions.
46 changes: 28 additions & 18 deletions WPS/geogrid/src/process_tile_module.F
Expand Up @@ -885,9 +885,9 @@ subroutine process_tile(which_domain, grid_type, dynopt, &
i1=field_count,i2=NUM_FIELDS-NUM_AUTOMATIC_FIELDS,s1=gradname)

if (grid_type == 'C') then
call calc_dfdx(field, slp_field, sm1, sm2, min_level, em1, em2, max_level, mapfac_ptr_x)
call calc_dfdx(field, slp_field, sm1, sm2, min_level, em1, em2, max_level, mapfac_ptr_x, sub_x)
else if (grid_type == 'E') then
call calc_dfdx(field, slp_field, sm1, sm2, min_level, em1, em2, max_level)
call calc_dfdx(field, slp_field, sm1, sm2, min_level, em1, em2, max_level, sr_x=sub_x)
end if
call write_field(sm1, em1, sm2, em2, &
min_level, max_level, trim(gradname), datestr, real_array=slp_field)
Expand All @@ -902,9 +902,9 @@ subroutine process_tile(which_domain, grid_type, dynopt, &
i1=field_count,i2=NUM_FIELDS-NUM_AUTOMATIC_FIELDS,s1=gradname)

if (grid_type == 'C') then
call calc_dfdy(field, slp_field, sm1, sm2, min_level, em1, em2, max_level, mapfac_ptr_y)
call calc_dfdy(field, slp_field, sm1, sm2, min_level, em1, em2, max_level, mapfac_ptr_y, sub_y)
else if (grid_type == 'E') then
call calc_dfdy(field, slp_field, sm1, sm2, min_level, em1, em2, max_level)
call calc_dfdy(field, slp_field, sm1, sm2, min_level, em1, em2, max_level, sr_y=sub_y)
end if
call write_field(sm1, em1, sm2, em2, &
min_level, max_level, trim(gradname), datestr, real_array=slp_field)
Expand Down Expand Up @@ -1942,7 +1942,7 @@ end subroutine process_neighbor
! the result in dst_array.
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine calc_dfdy(src_arr, dst_arr, start_mem_i, start_mem_j, start_mem_k, &
end_mem_i, end_mem_j, end_mem_k, mapfac)
end_mem_i, end_mem_j, end_mem_k, mapfac, sr_y)

! Modules
use gridinfo_module
Expand All @@ -1954,40 +1954,45 @@ subroutine calc_dfdy(src_arr, dst_arr, start_mem_i, start_mem_j, start_mem_k, &
real, dimension(start_mem_i:end_mem_i, start_mem_j:end_mem_j,start_mem_k:end_mem_k), intent(in) :: src_arr
real, dimension(start_mem_i:end_mem_i, start_mem_j:end_mem_j,start_mem_k:end_mem_k), intent(out) :: dst_arr
real, dimension(start_mem_i:end_mem_i, start_mem_j:end_mem_j), intent(in), optional :: mapfac
integer, intent(in), optional :: sr_y

! Local variables
integer :: i, j, k
real :: l_sr_y

l_sr_y=1.
if (present(sr_y)) l_sr_y=sr_y

if (present(mapfac)) then
do k=start_mem_k,end_mem_k
do i=start_mem_i, end_mem_i
do j=start_mem_j+1, end_mem_j-1
dst_arr(i,j,k) = (src_arr(i,j+1,k) - src_arr(i,j-1,k))/(2.*dykm*mapfac(i,j))
dst_arr(i,j,k) = sr_y*(src_arr(i,j+1,k) - src_arr(i,j-1,k))/(2.*dykm*mapfac(i,j))
end do
end do

do i=start_mem_i, end_mem_i
dst_arr(i,start_mem_j,k) = (src_arr(i,start_mem_j+1,k) - src_arr(i,start_mem_j,k))/(dykm*mapfac(i,j))
dst_arr(i,start_mem_j,k) = sr_y*(src_arr(i,start_mem_j+1,k) - src_arr(i,start_mem_j,k))/(dykm*mapfac(i,j))
end do

do i=start_mem_i, end_mem_i
dst_arr(i,end_mem_j,k) = (src_arr(i,end_mem_j,k) - src_arr(i,end_mem_j-1,k))/(dykm*mapfac(i,j))
dst_arr(i,end_mem_j,k) = sr_y*(src_arr(i,end_mem_j,k) - src_arr(i,end_mem_j-1,k))/(dykm*mapfac(i,j))
end do
end do
else
do k=start_mem_k,end_mem_k
do i=start_mem_i, end_mem_i
do j=start_mem_j+1, end_mem_j-1
dst_arr(i,j,k) = (src_arr(i,j+1,k) - src_arr(i,j-1,k))/(2.*dykm)
dst_arr(i,j,k) = sr_y*(src_arr(i,j+1,k) - src_arr(i,j-1,k))/(2.*dykm)
end do
end do

do i=start_mem_i, end_mem_i
dst_arr(i,start_mem_j,k) = (src_arr(i,start_mem_j+1,k) - src_arr(i,start_mem_j,k))/(dykm)
dst_arr(i,start_mem_j,k) = sr_y*(src_arr(i,start_mem_j+1,k) - src_arr(i,start_mem_j,k))/(dykm)
end do

do i=start_mem_i, end_mem_i
dst_arr(i,end_mem_j,k) = (src_arr(i,end_mem_j,k) - src_arr(i,end_mem_j-1,k))/(dykm)
dst_arr(i,end_mem_j,k) = sr_y*(src_arr(i,end_mem_j,k) - src_arr(i,end_mem_j-1,k))/(dykm)
end do
end do
end if
Expand All @@ -2002,7 +2007,7 @@ end subroutine calc_dfdy
! the result in dst_array.
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine calc_dfdx(src_arr, dst_arr, start_mem_i, start_mem_j, &
start_mem_k, end_mem_i, end_mem_j, end_mem_k, mapfac)
start_mem_k, end_mem_i, end_mem_j, end_mem_k, mapfac, sr_x)

! Modules
use gridinfo_module
Expand All @@ -2014,40 +2019,45 @@ subroutine calc_dfdx(src_arr, dst_arr, start_mem_i, start_mem_j, &
real, dimension(start_mem_i:end_mem_i, start_mem_j:end_mem_j, start_mem_k:end_mem_k), intent(in) :: src_arr
real, dimension(start_mem_i:end_mem_i, start_mem_j:end_mem_j, start_mem_k:end_mem_k), intent(out) :: dst_arr
real, dimension(start_mem_i:end_mem_i, start_mem_j:end_mem_j), intent(in), optional :: mapfac
integer, intent(in), optional :: sr_x

! Local variables
integer :: i, j, k
real :: l_sr_x

l_sr_x=1.
if (present(sr_x)) l_sr_x=sr_x

if (present(mapfac)) then
do k=start_mem_k, end_mem_k
do i=start_mem_i+1, end_mem_i-1
do j=start_mem_j, end_mem_j
dst_arr(i,j,k) = (src_arr(i+1,j,k) - src_arr(i-1,j,k))/(2.*dxkm*mapfac(i,j))
dst_arr(i,j,k) = sr_x*(src_arr(i+1,j,k) - src_arr(i-1,j,k))/(2.*dxkm*mapfac(i,j))
end do
end do

do j=start_mem_j, end_mem_j
dst_arr(start_mem_i,j,k) = (src_arr(start_mem_i+1,j,k) - src_arr(start_mem_i,j,k))/(dxkm*mapfac(i,j))
dst_arr(start_mem_i,j,k) = sr_x*(src_arr(start_mem_i+1,j,k) - src_arr(start_mem_i,j,k))/(dxkm*mapfac(i,j))
end do

do j=start_mem_j, end_mem_j
dst_arr(end_mem_i,j,k) = (src_arr(end_mem_i,j,k) - src_arr(end_mem_i-1,j,k))/(dxkm*mapfac(i,j))
dst_arr(end_mem_i,j,k) = sr_x*(src_arr(end_mem_i,j,k) - src_arr(end_mem_i-1,j,k))/(dxkm*mapfac(i,j))
end do
end do
else
do k=start_mem_k, end_mem_k
do i=start_mem_i+1, end_mem_i-1
do j=start_mem_j, end_mem_j
dst_arr(i,j,k) = (src_arr(i+1,j,k) - src_arr(i-1,j,k))/(2.*dxkm)
dst_arr(i,j,k) = sr_x*(src_arr(i+1,j,k) - src_arr(i-1,j,k))/(2.*dxkm)
end do
end do

do j=start_mem_j, end_mem_j
dst_arr(start_mem_i,j,k) = (src_arr(start_mem_i+1,j,k) - src_arr(start_mem_i,j,k))/(dxkm)
dst_arr(start_mem_i,j,k) = sr_x*(src_arr(start_mem_i+1,j,k) - src_arr(start_mem_i,j,k))/(dxkm)
end do

do j=start_mem_j, end_mem_j
dst_arr(end_mem_i,j,k) = (src_arr(end_mem_i,j,k) - src_arr(end_mem_i-1,j,k))/(dxkm)
dst_arr(end_mem_i,j,k) = sr_x*(src_arr(end_mem_i,j,k) - src_arr(end_mem_i-1,j,k))/(dxkm)
end do
end do
end if
Expand Down

0 comments on commit 80d70fb

Please sign in to comment.