da_tv_profile_adj.inc
References to this file elsewhere.
1 subroutine da_tv_profile_adj(xp,jo_grad_x,xb,i,j,dx,dxm,dy,dym, &
2 pre_ma,tv_ma,ADJ_pre_ma,ADJ_tv_ma)
3
4 !-----------------------------------------------------------------------
5 ! Purpose: adjoint routine for da_tv_profile
6 !-----------------------------------------------------------------------
7
8 implicit none
9
10 type (x_type) , intent(inout) :: jo_grad_x ! grad_x(jo)
11 type (xb_type), intent(in) :: xb ! first guess state.
12 type (xpose_type), intent(in) :: xp ! Dimensions and xpose buffers.
13 integer, intent(in) :: i, j ! OBS location
14 real, intent(in) :: dx, dxm ! interpolation weights.
15 real, intent(in) :: dy, dym ! interpolation weights.
16
17 real, dimension(xp%kts-1:xp%kte+1), intent(in) :: pre_ma,tv_ma
18 real, dimension(xp%kts-1:xp%kte+1), intent(inout) :: ADJ_pre_ma,ADJ_tv_ma
19
20 integer :: ii,jj,ks,ke
21 real, dimension(2,2,xp%kts:xp%kte) :: ADJ_tv_m
22
23 ks = xp%kts; ke = xp%kte
24
25 ADJ_tv_m(1,1,ks:ke) = dym*dxm * ADJ_tv_ma (ks:ke)
26 ADJ_tv_m(2,1,ks:ke) = dym*dx * ADJ_tv_ma (ks:ke)
27 ADJ_tv_m(1,2,ks:ke) = dy*dxm* ADJ_tv_ma (ks:ke)
28 ADJ_tv_m(2,2,ks:ke) = dy*dx* ADJ_tv_ma (ks:ke)
29 jo_grad_x%p(i,j,ks:ke) = jo_grad_x%p(i,j,ks:ke) + dym*dxm * ADJ_pre_ma(ks:ke)
30 jo_grad_x%p(i+1,j,ks:ke) = jo_grad_x%p(i+1,j,ks:ke) + dym*dx * ADJ_pre_ma(ks:ke)
31 jo_grad_x%p(i,j+1,ks:ke) = jo_grad_x%p(i,j+1,ks:ke) + dy*dxm * ADJ_pre_ma(ks:ke)
32 jo_grad_x%p(i+1,j+1,ks:ke)= jo_grad_x%p(i+1,j+1,ks:ke) + dy*dx* ADJ_pre_ma(ks:ke)
33 ADJ_tv_ma (ks:ke) = 0.
34 ADJ_pre_ma(ks:ke) = 0.
35
36 do ii=i,i+1
37 do jj=j,j+1
38 jo_grad_x%t(ii,jj,ks:ke) = jo_grad_x%t(ii,jj,ks:ke) + &
39 ADJ_tv_m(ii-i+1,jj-j+1,ks:ke)*(1.+0.61*xb%q(ii,jj,ks:ke))
40 jo_grad_x%q(ii,jj,ks:ke) = jo_grad_x%q(ii,jj,ks:ke) + &
41 0.61*xb%t(ii,jj,ks:ke)*ADJ_tv_m(ii-i+1,jj-j+1,ks:ke)
42 end do
43 end do
44
45 end subroutine da_tv_profile_adj
46
47