module omegas_ad_mod
!$$$ module documentation block
! . . . .
! module: omegas_ad_mod module wrapper around subroutine omegas_ad
! prgmmr: parrish org: np22 date: 2013-01-26
!
! abstract: This module has been added as a wrapper around subroutine omegas_ad
! to eliminate type mismatch compile errors when using the debug
! compile option on WCOSS.
!
! program history log:
! 2012-01-26 parrish
!
! subroutines included:
! omegas_ad
! omegas_ad_1_1_
! omegas_ad_im_ix_
implicit none
! set default to private
private
! set subroutines to public
public :: omegas_ad
interface omegas_ad
module procedure omegas_ad_1_1_
module procedure omegas_ad_im_ix_
end interface
contains
subroutine omegas_ad_1_1_( im, ix, km, dphi_i_, dlam_i_, u_i_, v_i_, div_i_, &
ps_i_, rcl_, del_, sl_, vvel_o_, u_i_ad_, v_i_ad_, div_i_ad_, vvel_o_ad_, &
adjoint)
!$$$ subprogram documentation block
! . . . .
! subprogram: omegas_ad_1_1_
! prgmmr: parrish org: np22 date: 2013-01-26
!
! abstract: interface for omegas_ad, where im=1,ix=1, and calling routine has
! no dimension index corresponding to im and ix.
!
! program history log:
! 2013-01-26 parrish - initial documentation
!
! input argument list:
! im,ix,km,adjoint
! dphi_i_,dlam_i_,u_i_(km),v_i_(km),div_i_(km),ps_i_,rcl_,del_(km),sl_(km)
! vvel_o_(km),u_i_ad_(km),v_i_ad_(km),div_i_ad_(km),vvel_o_ad_(km)
!
! output argument list:
! vvel_o_(km),u_i_ad_(km),v_i_ad_(km),div_i_ad_(km),vvel_o_ad_(km)
use kinds, only: r_kind,i_kind
implicit none
integer(i_kind),intent(in ) :: im,ix,km
logical ,intent(in ) :: adjoint
real(r_kind) ,intent(in ) :: dphi_i_,dlam_i_,u_i_(km),v_i_(km),div_i_(km)
real(r_kind) ,intent(in ) :: ps_i_,rcl_
real(r_kind) ,intent(in ) :: del_(km),sl_(km)
real(r_kind) ,intent(inout) :: vvel_o_(km),u_i_ad_(km),v_i_ad_(km)
real(r_kind) ,intent(inout) :: div_i_ad_(km),vvel_o_ad_(km)
real(r_kind) :: dphi_i(ix),dlam_i(ix),u_i(km,ix),v_i(km,ix),div_i(km,ix)
real(r_kind) :: ps_i(ix),rcl(ix)
real(r_kind) :: del(km,ix),sl(km,ix)
real(r_kind) :: vvel_o(km,ix),u_i_ad(km,ix),v_i_ad(km,ix)
real(r_kind) :: div_i_ad(km,ix),vvel_o_ad(km,ix)
integer(i_kind) k
if( im /= 1 .or. ix /= 1 ) then
write(6,*)' GSCOND_AD_1_1_, IM,IX=',IM,IX,' -- BOTH MUST BE 1. PROGRAM FAILS'
stop
end if
dphi_i(1)=dphi_i_
dlam_i(1)=dlam_i_
ps_i(1)=ps_i_
rcl(1)=rcl_
do k=1,km
u_i(k,1)=u_i_(k)
v_i(k,1)=v_i_(k)
div_i(k,1)=div_i_(k)
del(k,1)=del_(k)
sl(k,1)=sl_(k)
vvel_o(k,1)=vvel_o_(k)
u_i_ad(k,1)=u_i_ad_(k)
v_i_ad(k,1)=v_i_ad_(k)
div_i_ad(k,1)=div_i_ad_(k)
vvel_o_ad(k,1)=vvel_o_ad_(k)
end do
call omegas_ad_im_ix_( im, ix, km, dphi_i, dlam_i, u_i, v_i, div_i, &
ps_i, rcl, del, sl, vvel_o, u_i_ad, v_i_ad, div_i_ad, vvel_o_ad, &
adjoint)
do k=1,km
vvel_o_(k)=vvel_o(k,1)
u_i_ad_(k)=u_i_ad(k,1)
v_i_ad_(k)=v_i_ad(k,1)
div_i_ad_(k)=div_i_ad(k,1)
vvel_o_ad_(k)=vvel_o_ad(k,1)
end do
end subroutine omegas_ad_1_1_
subroutine omegas_ad_im_ix_( im, ix, km, dphi_i, dlam_i, u_i, v_i, div_i, &
ps_i, rcl, del, sl, vvel_o, u_i_ad, v_i_ad, div_i_ad, vvel_o_ad, &
adjoint)
!$$$ subprogram documentation block
! . . . .
! subprogram: omegas_ad forward & adjoint model for GFS vertical velocity calculation
! prgmmr: treadon org: np23 date: 2003-12-18
!
! abstract: This subroutine contains the forward and ajoint models for the
! GFS vertical velocity calculation
!
! program history log:
! 2003-12-18 treadon - initial routine
! 2004-06-15 treadon - reformat documenation
! 2006-04-12 treadon - change del and sl from 1d to 2d arrays
! 2006-08-10 treadon - change dphi_i,dlam_i from ln(ps) to ps
! 2013-01-26 parrish - module added as a wrapper around subroutine omegas_ad
! to eliminate type mismatch compile errors when
! using the debug
! compile option on WCOSS.
!
! input argument list:
! im - actual number of profiles to be processed
! ix - maximum number of profiles to process (array dimension)
! km - number of levels in vertical profile
! dphi_i - partial derivative of ps with respect to latitude
! dlam_i - partial derivative of ps with respect to longitude
! u_i - zonal wind
! v_i - meridional wind
! div_i - horizontal divergence
! ps_i - surface pressure
! rcl - 1/sin(lat)**2
! del - "sigma" thickness of layers
! sl - "sigma" value at layer midpoints
! vvel_o_ad- vertical velocity perturbation
! adjoint - logical flag (.false.=forward model only, .true.=forward and ajoint)
!
! output argument list:
! vvel_o - vertical velocity
! vvel_o_ad- vertical velocity perturbation
! u_i_ad - partial derivative of vertical velocity with respect to zonal wind
! v_i_ad - partial derivative of vertical velocity with respect to meridional wind
! div_i_ad - partial derivative of vertical velocity with respect to horizontal divergence
!
!
! remarks:
! The adjoint portion of this routine was generated by the
! Tangent linear and Adjoint Model Compiler, TAMC 5.3.0
!
! attributes:
! language: f90
! machine: ibm RS/6000 SP
!
!$$$
!==============================================
! all entries are defined explicitly
!==============================================
use kinds, only: r_kind,i_kind
use constants, only: zero, half
implicit none
!==============================================
! define arguments
!==============================================
logical ,intent(in ) :: adjoint
integer(i_kind),intent(in ) :: ix
integer(i_kind),intent(in ) :: km
real(r_kind) ,intent( out) :: div_i_ad(km,ix)
real(r_kind) ,intent( out) :: u_i_ad(km,ix)
real(r_kind) ,intent( out) :: v_i_ad(km,ix)
real(r_kind) ,intent(inout) :: vvel_o_ad(km,ix)
real(r_kind) ,intent(in ) :: del(km,ix)
real(r_kind) ,intent(in ) :: div_i(km,ix)
real(r_kind) ,intent(in ) :: dlam_i(ix)
real(r_kind) ,intent(in ) :: dphi_i(ix)
integer(i_kind),intent(in ) :: im
real(r_kind) ,intent(in ) :: ps_i(ix)
real(r_kind) ,intent(in ) :: rcl(ix)
real(r_kind) ,intent(in ) :: sl(km,ix)
real(r_kind) ,intent(in ) :: u_i(km,ix)
real(r_kind) ,intent(in ) :: v_i(km,ix)
real(r_kind) ,intent( out) :: vvel_o(km,ix)
!==============================================
! define local variables
!==============================================
real(r_kind) cb_ad(km,ix)
real(r_kind) cg_ad(km,ix)
real(r_kind) db_ad(km,ix)
real(r_kind) dot_ad(km+1,ix)
real(r_kind) cb(km,ix)
real(r_kind) cg(km,ix)
real(r_kind) db(km,ix)
real(r_kind) dlam(ix)
real(r_kind) dot(km+1,ix)
real(r_kind) dphi(ix)
integer(i_kind) i
integer(i_kind) ip1
integer(i_kind) ip2
integer(i_kind) k
!----------------------------------------------
! RESET LOCAL ADJOINT VARIABLES
!----------------------------------------------
do ip2 = 1, ix
do ip1 = 1, km
cb_ad(ip1,ip2) = zero
end do
end do
do ip2 = 1, ix
do ip1 = 1, km
cg_ad(ip1,ip2) = zero
end do
end do
do ip2 = 1, ix
do ip1 = 1, km
db_ad(ip1,ip2) = zero
end do
end do
do ip2 = 1, ix
do ip1 = 1, km+1
dot_ad(ip1,ip2) = zero
end do
end do
!----------------------------------------------
! ROUTINE BODY
!----------------------------------------------
!----------------------------------------------
! FUNCTION AND TAPE COMPUTATIONS
!----------------------------------------------
do i = 1, im
do k = 1, km+1
dot(k,i) = zero
end do
end do
do i = 1, im
dphi(i) = dphi_i(i)*rcl(i) / ps_i(i)
dlam(i) = dlam_i(i)*rcl(i) / ps_i(i)
end do
do i = 1, im
do k = 1, km
cg(k,i) = u_i(k,i)*dlam(i)+v_i(k,i)*dphi(i)
end do
end do
do i = 1, im
db(1,i) = del(1,i)*div_i(1,i)
cb(1,i) = del(1,i)*cg(1,i)
end do
do i = 1, im
do k = 1, km-1
db(k+1,i) = db(k,i)+del(k+1,i)*div_i(k+1,i)
cb(k+1,i) = cb(k,i)+del(k+1,i)*cg (k+1,i)
end do
end do
do i = 1, im
do k = 1, km-1
dot(k+1,i) = dot(k,i)+del(k,i)*(db(km,i)+cb(km,i)-div_i(k,i)- &
cg(k,i))
end do
end do
do i = 1, im
do k = 1, km
vvel_o(k,i) = sl(k,i)*(cg(k,i)-cb(km,i)-db(km,i))-half*(dot(k+1, &
i)+dot(k,i))
vvel_o(k,i) = vvel_o(k,i)*ps_i(i)
end do
end do
if (.not.adjoint) return
!----------------------------------------------
! ADJOINT COMPUTATIONS
!----------------------------------------------
do i = 1, im
do k = 1, km
vvel_o_ad(k,i) = vvel_o_ad(k,i)*ps_i(i)
cb_ad(km,i) = cb_ad(km,i)-vvel_o_ad(k,i)*sl(k,i)
cg_ad(k,i) = cg_ad(k,i)+vvel_o_ad(k,i)*sl(k,i)
db_ad(km,i) = db_ad(km,i)-vvel_o_ad(k,i)*sl(k,i)
dot_ad(k+1,i) = dot_ad(k+1,i)-half*vvel_o_ad(k,i)
dot_ad(k,i) = dot_ad(k,i)-half*vvel_o_ad(k,i)
vvel_o_ad(k,i) = zero
end do
end do
do i = 1, im
do k = km-1, 1, -1
cb_ad(km,i) = cb_ad(km,i)+dot_ad(k+1,i)*del(k,i)
cg_ad(k,i) = cg_ad(k,i)-dot_ad(k+1,i)*del(k,i)
db_ad(km,i) = db_ad(km,i)+dot_ad(k+1,i)*del(k,i)
div_i_ad(k,i) = div_i_ad(k,i)-dot_ad(k+1,i)*del(k,i)
dot_ad(k,i) = dot_ad(k,i)+dot_ad(k+1,i)
dot_ad(k+1,i) = zero
end do
end do
do i = 1, im
do k = km-1, 1, -1
cg_ad(k+1,i) = cg_ad(k+1,i)+cb_ad(k+1,i)*del(k+1,i)
cb_ad(k,i) = cb_ad(k,i)+cb_ad(k+1,i)
cb_ad(k+1,i) = zero
div_i_ad(k+1,i) = div_i_ad(k+1,i)+db_ad(k+1,i)*del(k+1,i)
db_ad(k,i) = db_ad(k,i)+db_ad(k+1,i)
db_ad(k+1,i) = zero
end do
end do
do i = 1, im
cg_ad(1,i) = cg_ad(1,i)+cb_ad(1,i)*del(1,i)
cb_ad(1,i) = zero
div_i_ad(1,i) = div_i_ad(1,i)+db_ad(1,i)*del(1,i)
db_ad(1,i) = zero
end do
do i = 1, im
do k = 1, km
u_i_ad(k,i) = u_i_ad(k,i)+cg_ad(k,i)*dlam(i)
v_i_ad(k,i) = v_i_ad(k,i)+cg_ad(k,i)*dphi(i)
cg_ad(k,i) = zero
end do
end do
return
end subroutine omegas_ad_im_ix_
end module omegas_ad_mod