da_transfer_xatowrftl_adj.inc
References to this file elsewhere.
1 subroutine da_transfer_xatowrftl_adj(grid, config_flags, filnam)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Convert WRFTL variables to analysis increments
5 ! (inverse of the incremental part of xatowrf)
6 !---------------------------------------------------------------------------
7
8 implicit none
9
10 type(domain), intent(inout) :: grid
11 type(grid_config_rec_type), intent(inout) :: config_flags
12
13 character*4, intent(in) :: filnam
14
15 ! Local variables
16
17 integer :: i, j, k
18 integer :: is, ie, js, je, ks, ke
19 real :: sdmd, s1md
20 real :: rho_cgrid
21
22 real, dimension(grid%xp%ims:grid%xp%ime,grid%xp%jms:grid%xp%jme,grid%xp%kms:grid%xp%kme) :: a_press
23
24 integer ndynopt
25
26 if (trace_use) call da_trace_entry("da_transfer_xatowrftl_adj")
27
28 is=grid%xp%its
29 ie=grid%xp%ite
30 js=grid%xp%jts
31 je=grid%xp%jte
32 ks=grid%xp%kts
33 ke=grid%xp%kte
34
35 !---------------------------------------------------------------------------
36 ! [7.0] Adjoint of outPUT (inPUT)
37 !---------------------------------------------------------------------------
38
39 ndynopt = grid%dyn_opt
40 grid%dyn_opt = DYN_EM_AD
41 call nl_set_dyn_opt (1 , DYN_EM_AD)
42
43 call da_med_initialdata_input(grid , config_flags, filnam)
44
45 grid%dyn_opt = ndynopt
46 call nl_set_dyn_opt (1 , DYN_EM)
47
48 !---------------------------------------------------------------------------
49 ! [6.0] Adjoint of save OTHERinCREMENT
50 !---------------------------------------------------------------------------
51
52 do k=ks,ke+1
53 do j=js,je
54 do i=is,ie
55 grid%xa%w(i,j,k)=grid%em_a_w_2(i,j,k)
56 end do
57 end do
58 end do
59
60 #ifdef VAR4D_MICROPHYSICS
61 ! New code needed once we introduce the microphysics to 4dvar in 2008
62 if (size(grid%moist,dim=4) >= 4) then
63 do k=ks,ke
64 do j=js,je
65 do i=is,ie
66 grid%xa%qcw(i,j,k)=grid%a_moist(i,j,k,p_qcw)
67 grid%xa%qrn(i,j,k)=grid%a_moist(i,j,k,p_qrn)
68 end do
69 end do
70 end do
71 end if
72
73 if (size(grid%moist,dim=4) >= 6) then
74 do k=ks,ke
75 do j=js,je
76 do i=is,ie
77 grid%xa%qci(i,j,k)=grid%a_moist(i,j,k,p_qci)
78 grid%xa%qsn(i,j,k)=grid%a_moist(i,j,k,p_qsn)
79 end do
80 end do
81 end do
82 end if
83
84 if (size(grid%moist,dim=4) >= 7) then
85 do k=ks,ke
86 do j=js,je
87 do i=is,ie
88 grid%xa%qgr(i,j,k)=grid%a_moist(i,j,k,p_qgr)
89 end do
90 end do
91 end do
92 end if
93
94 #endif
95
96
97 !---------------------------------------------------------------------------
98 ! [5.0] Adjoint of CONVERT FROM A-GRID TO C-GRID
99 !---------------------------------------------------------------------------
100
101 ! Fill the halo region for u and v.
102 call wrf_dm_halo(grid%xp%domdesc,grid%xp%comms,grid%xp%halo_id3)
103
104 do k=ks,ke
105 do j=js,je
106 do i=is,ie
107 grid%xa%u(i,j,k)=0.5*(grid%em_a_u_2(i+1,j ,k)+grid%em_a_u_2(i,j,k))
108 grid%xa%v(i,j,k)=0.5*(grid%em_a_v_2(i ,j+1,k)+grid%em_a_v_2(i,j,k))
109 end do
110 end do
111 end do
112
113 #ifdef DM_PARALLEL
114
115 ! The western boundary
116 if (is == grid%xp%ids) then
117 grid%xa%u(is ,js:je,ks:ke)=grid%xa%u(is ,js:je,ks:ke)+grid%em_a_u_2(is,js:je,ks:ke)
118 grid%xa%u(is+1,js:je,ks:ke)=grid%xa%u(is+1,js:je,ks:ke)-grid%em_a_u_2(is,js:je,ks:ke)/2.
119 end if
120
121 ! The eastern boundary
122 if (ie == grid%xp%ide) then
123 grid%xa%u(ie ,js:je,ks:ke)=grid%xa%u(ie ,js:je,ks:ke)+grid%em_a_u_2(ie+1,js:je,ks:ke)
124 grid%xa%u(ie-1,js:je,ks:ke)=grid%xa%u(ie-1,js:je,ks:ke)-grid%em_a_u_2(ie+1,js:je,ks:ke)/2.
125 end if
126
127 ! The southern boundary
128 if (js == grid%xp%jds) then
129 grid%xa%v(is:ie,js ,ks:ke)=grid%xa%v(is:ie,js ,ks:ke)+grid%em_a_v_2(is:ie,js,ks:ke)
130 grid%xa%v(is:ie,js+1,ks:ke)=grid%xa%v(is:ie,js+1,ks:ke)-grid%em_a_v_2(is:ie,js,ks:ke)/2.
131 end if
132
133 ! The northern boundary
134 if (je == grid%xp%jde) then
135 grid%xa%v(is:ie,je ,ks:ke)=grid%xa%v(is:ie,je ,ks:ke)+grid%em_a_v_2(is:ie,je+1,ks:ke)
136 grid%xa%v(is:ie,js-1,ks:ke)=grid%xa%v(is:ie,je-1,ks:ke)-grid%em_a_v_2(is:ie,je+1,ks:ke)/2.
137 end if
138
139 #else
140
141 !The western boundary
142 grid%xa%u(is ,js:je,ks:ke)=grid%xa%u(is ,js:je,ks:ke)+grid%em_a_u_2(is,js:je,ks:ke)
143 grid%xa%u(is+1,js:je,ks:ke)=grid%xa%u(is+1,js:je,ks:ke)-grid%em_a_u_2(is,js:je,ks:ke)/2.
144
145 ! The eastern boundary
146 grid%xa%u(ie ,js:je,ks:ke)=grid%xa%u(ie ,js:je,ks:ke)+grid%em_a_u_2(ie+1,js:je,ks:ke)
147 grid%xa%u(ie-1,js:je,ks:ke)=grid%xa%u(ie-1,js:je,ks:ke)-grid%em_a_u_2(ie+1,js:je,ks:ke)/2.
148
149 ! The southern boundary
150 grid%xa%v(is:ie,js ,ks:ke)=grid%xa%v(is:ie,js ,ks:ke)+grid%em_a_v_2(is:ie,js,ks:ke)
151 grid%xa%v(is:ie,js+1,ks:ke)=grid%xa%v(is:ie,js+1,ks:ke)-grid%em_a_v_2(is:ie,js,ks:ke)/2.
152
153 ! The northern boundary
154 grid%xa%v(is:ie,je ,ks:ke)=grid%xa%v(is:ie,je ,ks:ke)+grid%em_a_v_2(is:ie,je+1,ks:ke)
155 grid%xa%v(is:ie,js-1,ks:ke)=grid%xa%v(is:ie,je-1,ks:ke)-grid%em_a_v_2(is:ie,je+1,ks:ke)/2.
156
157 #endif
158
159 !---------------------------------------------------------------------------
160 ! [4.0] Adjoint of CONVERT TEMPERATURE inCREMENTS inTO THETA inCREMENTS
161 ! EVALUATE ALSO THE inCREMENTS OF (1/rho) AND GEOPOTENTIAL
162 !---------------------------------------------------------------------------
163
164 a_press(is:ie,js:je,ks:ke+1)=0.
165 do k=ke,ks,-1
166 do j=js,je
167 do i=is,ie
168 rho_cgrid=-(grid%em_ph_2(i,j,k+1)-grid%em_ph_2(i,j,k))*grid%em_a_ph_2(i,j,k+1)/grid%xb%rho(i,j,k)
169 a_press(i,j,k )=a_press(i,j,k )+grid%em_a_ph_2(i,j,k+1)/grid%xb%rho(i,j,k)
170 a_press(i,j,k+1)=a_press(i,j,k+1)-grid%em_a_ph_2(i,j,k+1)/grid%xb%rho(i,j,k)
171 grid%em_a_ph_2(i,j,k ) =grid%em_a_ph_2(i,j,k) +grid%em_a_ph_2(i,j,k+1)
172 grid%xa%q(i,j,k)=-grid%xb%rho(i,j,k)*0.61*rho_cgrid/(1.+0.61*grid%xb%q(i,j,k))
173 grid%xa%t(i,j,k)=-grid%xb%rho(i,j,k)*rho_cgrid/grid%xb%t(i,j,k)
174 grid%xa%p(i,j,k)= grid%xb%rho(i,j,k)*rho_cgrid/grid%xb%p(i,j,k)
175 end do
176 end do
177 end do
178
179 do k=ks,ke
180 do j=js,je
181 do i=is,ie
182 grid%xa%p(i,j,k)=grid%xa%p(i,j,k)-(300.0+grid%em_t_2(i,j,k))*kappa*grid%em_a_t_2(i,j,k)/grid%xb%p(i,j,k)
183 grid%xa%t(i,j,k)=grid%xa%t(i,j,k)+(300.0+grid%em_t_2(i,j,k))*grid%em_a_t_2(i,j,k)/grid%xb%t(i,j,k)
184 end do
185 end do
186 end do
187
188 !---------------------------------------------------------------------------
189 ! [3.0] Adjoint of COMPUTE pressure increments (for computing theta increments)
190 !---------------------------------------------------------------------------
191
192 do k=ks,ke
193 do j=js,je
194 do i=is,ie
195 a_press(i,j,k+1)=a_press(i,j,k+1)+0.5*grid%xa%p(i,j,k)
196 a_press(i,j,k )=a_press(i,j,k )+0.5*grid%xa%p(i,j,k)
197 grid%xa%p(i,j,k)=0.
198 grid%a_moist(i,j,k,P_A_QV)=grid%a_moist(i,j,k,P_A_QV)-(grid%em_mu_2(i,j)+grid%em_mub(i,j))*a_press(i,j,k)*grid%em_dn(k)
199 grid%em_a_mu_2(i,j)=grid%em_a_mu_2(i,j)-a_press(i,j,k)*(1.0+grid%moist(i,j,k,P_QV))*grid%em_dn(k)
200 a_press(i,j,k+1)=a_press(i,j,k+1)+a_press(i,j,k)
201 end do
202 end do
203 end do
204
205 !---------------------------------------------------------------------------
206 ! [2.0] Adjoint of COMPUTE increments of dry-column air mass per unit area
207 !---------------------------------------------------------------------------
208
209 do j=js,je
210 do i=is,ie
211 sdmd=0.0
212 s1md=0.0
213 do k=ks,ke
214 s1md=s1md+(1.0+grid%moist(i,j,k,P_QV))*grid%em_dnw(k)
215 end do
216 sdmd=sdmd-grid%xb%psac(i,j)*grid%em_a_mu_2(i,j)/s1md
217 grid%xa%psfc(i,j)=-grid%em_a_mu_2(i,j)/s1md
218 do k=ks,ke
219 grid%a_moist(i,j,k,P_A_QV)=grid%a_moist(i,j,k,P_A_QV)+sdmd*grid%em_dnw(k)
220 end do
221 end do
222 end do
223
224 !---------------------------------------------------------------------------
225 ! [1.0] Adjoint of Get the mixing ratio of moisture
226 !---------------------------------------------------------------------------
227
228 do k=ks,ke
229 do j=js,je
230 do i=is,ie
231 grid%xa%q(i,j,k)=grid%xa%q(i,j,k)+grid%a_moist(i,j,k,P_A_QV)/(1.0-grid%xb%q(i,j,k))**2
232 end do
233 end do
234 end do
235
236 if (trace_use) call da_trace_exit("da_transfer_xatowrftl_adj")
237
238 end subroutine da_transfer_xatowrftl_adj
239
240