FV3 Bundle
wind_variables_mod.f90
Go to the documentation of this file.
1 ! (C) Copyright 2018 UCAR
2 !
3 ! This software is licensed under the terms of the Apache Licence Version 2.0
4 ! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
5 
6 !> Variable transforms on wind variables for fv3-jedi
7 !> Daniel Holdaway, NASA/JCSDA
8 
9 module wind_vt_mod
10 
11 use fv3jedi_constants_mod, only: pi, rad2deg
12 use fv3jedi_geom_mod, only: fv3jedi_geom
13 use fv3jedi_kinds_mod, only: kind_real
14 
15 implicit none
16 private
17 
18 public sfc_10m_winds
19 public uv_to_vortdivg
20 public vortdivg_to_psichi
21 public gauss_seidel
22 public psichi_to_uava
23 public psichi_to_uava_adm
24 public psichi_to_udvd
25 public d_to_ac
26 public d2a2c_vect
27 public divergence_corner
28 public a2b_ord4
29 public psichi_to_udvd_adm
30 public a2b_ord4_adm
31 public extrap_corner_adm
32 public a2d
33 public a2d_ad
34 public d2a
35 public d2a_ad
36 
37 contains
38 
39 !----------------------------------------------------------------------------
40 ! Lowest model level winds to 10m speed and direction -----------------------
41 !----------------------------------------------------------------------------
42 
43 subroutine sfc_10m_winds(geom,usrf,vsrf,f10r,spd10m,dir10m)
44 
45  implicit none
46 
47  !Arguments
48  type(fv3jedi_geom) , intent(in ) :: geom !Geometry for the model
49  real(kind=kind_real), intent(in ) :: usrf(geom%isc:geom%iec,geom%jsc:geom%jec) !Lowest model level u m/s
50  real(kind=kind_real), intent(in ) :: vsrf(geom%isc:geom%iec,geom%jsc:geom%jec) !Lowest model level v m/s
51  real(kind=kind_real), intent(in ) :: f10r(geom%isc:geom%iec,geom%jsc:geom%jec) !Ratio of lowest level to 10m
52  real(kind=kind_real), intent(out) :: spd10m(geom%isc:geom%iec,geom%jsc:geom%jec) !10m wind speed u m/s
53  real(kind=kind_real), intent(out) :: dir10m(geom%isc:geom%iec,geom%jsc:geom%jec) !10m model wind direction
54 
55  !Locals
56  integer :: isc,iec,jsc,jec,i,j
57  integer :: iquad
58  real(kind=kind_real) :: windangle, windratio
59  real(kind=kind_real), parameter :: windscale = 999999.0_kind_real
60  real(kind=kind_real), parameter :: windlimit = 0.0001_kind_real
61  real(kind=kind_real), parameter :: quadcof(4,2) = reshape((/ 0.0_kind_real, 1.0_kind_real, 1.0_kind_real, 2.0_kind_real, &
62  1.0_kind_real, -1.0_kind_real, 1.0_kind_real, -1.0_kind_real /), (/4, 2/))
63 
64  !In GSI these calculations are done after interpolation to obs location
65 
66  isc = geom%isc
67  iec = geom%iec
68  jsc = geom%jsc
69  jec = geom%jec
70 
71  !10m wind speed
72  spd10m(isc:iec,jsc:jec) = f10r(isc:iec,jsc:jec)*sqrt( usrf(isc:iec,jsc:jec)*usrf(isc:iec,jsc:jec) + &
73  vsrf(isc:iec,jsc:jec)*vsrf(isc:iec,jsc:jec) )
74 
75  !10m wind direction
76  do j = jsc,jec
77  do i = isc,iec
78 
79  if (usrf(i,j)*f10r(i,j) >= 0.0_kind_real .and. vsrf(i,j)*f10r(i,j) >= 0.0_kind_real) iquad = 1
80  if (usrf(i,j)*f10r(i,j) >= 0.0_kind_real .and. vsrf(i,j)*f10r(i,j) < 0.0_kind_real) iquad = 2
81  if (usrf(i,j)*f10r(i,j) < 0.0_kind_real .and. vsrf(i,j)*f10r(i,j) >= 0.0_kind_real) iquad = 3
82  if (usrf(i,j)*f10r(i,j) < 0.0_kind_real .and. vsrf(i,j)*f10r(i,j) < 0.0_kind_real) iquad = 4
83 
84  if (abs(vsrf(i,j)*f10r(i,j)) >= windlimit) then
85  windratio = (usrf(i,j)*f10r(i,j)) / (vsrf(i,j)*f10r(i,j))
86  else
87  windratio = 0.0_kind_real
88  if (abs(usrf(i,j)*f10r(i,j)) > windlimit) then
89  windratio = windscale * usrf(i,j)*f10r(i,j)
90  endif
91  endif
92 
93  windangle = atan(abs(windratio))
94 
95  dir10m(i,j) = rad2deg*(quadcof(iquad,1) * pi + windangle * quadcof(iquad, 2))
96 
97  enddo
98  enddo
99 
100 end subroutine sfc_10m_winds
101 
102 !----------------------------------------------------------------------------
103 
104 subroutine uv_to_vortdivg(geom,u,v,ua,va,vort,divg)
105 
106  use mpp_domains_mod, only: mpp_update_domains, mpp_get_boundary, dgrid_ne, cgrid_ne
107 
108  implicit none
109  type(fv3jedi_geom), intent(inout) :: geom
110  real(kind=kind_real), intent(inout) :: u(geom%isd:geom%ied ,geom%jsd:geom%jed+1,1:geom%npz) !Dgrid winds (u)
111  real(kind=kind_real), intent(inout) :: v(geom%isd:geom%ied+1,geom%jsd:geom%jed ,1:geom%npz) !Dgrid winds (v)
112  real(kind=kind_real), intent(inout) :: vort(geom%isd:geom%ied ,geom%jsd:geom%jed ,1:geom%npz) !Vorticity
113  real(kind=kind_real), intent(inout) :: divg(geom%isd:geom%ied ,geom%jsd:geom%jed ,1:geom%npz) !Divergence
114 
115  real(kind=kind_real), intent(inout) :: ua(geom%isd:geom%ied ,geom%jsd:geom%jed ,1:geom%npz)
116  real(kind=kind_real), intent(inout) :: va(geom%isd:geom%ied ,geom%jsd:geom%jed ,1:geom%npz)
117  real(kind=kind_real) :: uc(geom%isd:geom%ied+1,geom%jsd:geom%jed ,1:geom%npz)
118  real(kind=kind_real) :: vc(geom%isd:geom%ied ,geom%jsd:geom%jed+1,1:geom%npz)
119  real(kind=kind_real) :: ut(geom%isd:geom%ied+1,geom%jsd:geom%jed )
120  real(kind=kind_real) :: vt(geom%isd:geom%ied ,geom%jsd:geom%jed+1)
121 
122  real(kind=kind_real) :: wbuffer(geom%jsc:geom%jec,geom%npz)
123  real(kind=kind_real) :: sbuffer(geom%isc:geom%iec,geom%npz)
124  real(kind=kind_real) :: ebuffer(geom%jsc:geom%jec,geom%npz)
125  real(kind=kind_real) :: nbuffer(geom%isc:geom%iec,geom%npz)
126 
127  integer :: i,j,k,isc,iec,jsc,jec,npz
128  real(kind=kind_real) :: dt
129 
130  ! x-----------------x
131  ! | |
132  ! | |
133  ! | |
134  ! vd,uc| x |
135  ! | vort,divg |
136  ! | |
137  ! | |
138  ! x-----------------x
139  ! ud,vc
140 
141  isc = geom%isc
142  iec = geom%iec
143  jsc = geom%jsc
144  jec = geom%jec
145  npz = geom%npz
146 
147  ua = 0.0_kind_real
148  va = 0.0_kind_real
149 
150 
151  !Fill the halos of Dgrid winds and get A/C grid winds
152  !----------------------------------------------------
153  call d_to_ac(geom,u,v,ua,va,uc,vc)
154 
155 
156  !Calculate vorticity (Agrid)
157  !---------------------------
158 ! vort = 0.0_kind_real
159 ! do k=1,npz
160 ! do j=jsc,jec
161 ! do i=isc,iec
162 ! vort(i,j,k) = geom%rarea(i,j)*(u(i,j,k)*geom%dx(i,j)-u(i,j+1,k)*geom%dx(i,j+1) - &
163 ! v(i,j,k)*geom%dy(i,j)+v(i+1,j,k)*geom%dy(i+1,j))
164 ! enddo
165 ! enddo
166 ! enddo
167 
168  vort = 0.0_kind_real
169  do k=1,npz
170  do j=jsc,jec
171  do i=isc,iec
172  vort(i,j,k) = (va(i+1,j,k) - va(i-1,j,k))/(geom%dxc(i,j)+geom%dxc(i-1,j)) - (ua(i,j+1,k) - ua(i,j-1,k))/(geom%dyc(i,j)+geom%dyc(i,j-1))
173  enddo
174  enddo
175  enddo
176 
177  !Calculate divergence (Agrid)
178  !----------------------------
179  !Using C grid winds to give divergence at the grid center
180 ! divg = 0.0_kind_real
181 ! do k=1,npz
182 ! do j=jsc,jec
183 ! do i=isc,iec
184 ! divg(i,j,k) = geom%rarea(i,j)*(uc(i,j,k)*geom%dy(i,j)-uc(i,j+1,k)*geom%dy(i,j+1) + &
185 ! vc(i,j,k)*geom%dx(i,j)+vc(i+1,j,k)*geom%dx(i+1,j))
186 ! enddo
187 ! enddo
188 ! enddo
189 
190  divg = 0.0_kind_real
191  do k=1,npz
192  do j=jsc,jec
193  do i=isc,iec
194  divg(i,j,k) = (ua(i+1,j,k) - ua(i-1,j,k))/(geom%dxc(i,j)+geom%dxc(i-1,j)) + (va(i,j+1,k) - va(i,j-1,k))/(geom%dyc(i,j)+geom%dyc(i,j-1))
195  enddo
196  enddo
197  enddo
198 
199 endsubroutine uv_to_vortdivg
200 
201 !----------------------------------------------------------------------------
202 
203 subroutine vortdivg_to_psichi(geom,vort,divg,psi,chi)
204 
205  use mpi
206  use mpp_domains_mod, only: mpp_update_domains
207 
208  implicit none
209  type(fv3jedi_geom), intent(inout) :: geom
210  real(kind=kind_real), intent(in) :: vort(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz) !Vorticity
211  real(kind=kind_real), intent(in) :: divg(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz) !Divergence
212  real(kind=kind_real), intent(out) :: psi (geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz) !Stream function
213  real(kind=kind_real), intent(out) :: chi (geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz) !Velocity potential
214 
215  real(kind=kind_real) :: psinew(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz) !Stream function
216 
217  integer :: i,j,k,isc,iec,jsc,jec,npz
218  integer :: maxiter, iter, ierr
219  real(kind=kind_real) :: tolerance, converged, convergedg
220 
221  ! x-----------------x
222  ! | |
223  ! | |
224  ! | |
225  ! | x |
226  ! | vort,psi |
227  ! | divg,chi |
228  ! | |
229  ! x-----------------x
230 
231  call gauss_seidel(geom,psi,-vort)
232  call gauss_seidel(geom,chi,-divg)
233 
234 endsubroutine vortdivg_to_psichi
235 
236 !----------------------------------------------------------------------------
237 
238 subroutine gauss_seidel(geom,x,b)
239 
240  use mpi
241  use mpp_domains_mod, only: mpp_update_domains
242 
243  implicit none
244  type(fv3jedi_geom), intent(inout) :: geom
245  real(kind=kind_real), intent(inout) :: x(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz)
246  real(kind=kind_real), intent(in ) :: b(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz)
247 
248  integer :: i,j,k,isc,iec,jsc,jec,npz
249  integer :: maxiter, iter, ierr
250  real(kind=kind_real) :: tolerance, converged, convergedg
251  real(kind=kind_real) :: xnew(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz)
252 
253  isc = geom%isc
254  iec = geom%iec
255  jsc = geom%jsc
256  jec = geom%jec
257  npz = geom%npz
258 
259  x = 0.0_kind_real
260 
261  converged = 1.0_kind_real
262  tolerance = 10e-4_kind_real
263 
264  maxiter = 10000
265  iter = 0
266 
267  do while (iter < maxiter .and. converged > tolerance)
268 
269  iter = iter + 1
270 
271  xnew(isc:iec,jsc:jec,:) = x(isc:iec,jsc:jec,:)
272 
273  call mpp_update_domains(xnew, geom%domain, complete=.true.)
274 
275  !Rearranged Laplacian (not the corret dx's)
276  do k=1,npz
277  do j=jsc,jec
278  do i=isc,iec
279  xnew(i,j,k) = 0.25_kind_real*(xnew(i+1,j,k)+xnew(i,j+1,k)+xnew(i-1,j,k)+xnew(i,j-1,k)-geom%dxa(i,j)**2*b(i,j,k))
280  enddo
281  enddo
282  enddo
283 
284  !Check for convergence
285  converged = maxval(abs((xnew-x)))
286  call mpi_allreduce(converged,convergedg,1,mpi_real8,mpi_max,mpi_comm_world,ierr)
287 
288  !Print info
289  !print*, 'Gauss-Seidel iteration number and max difference, ', iter, convergedg
290 
291  !Update guess
292  x = xnew
293 
294  enddo
295 
296 endsubroutine gauss_seidel
297 
298 !----------------------------------------------------------------------------
299 
300 subroutine psichi_to_uava(geom,psi,chi,ua,va)
301 
302  use mpp_domains_mod, only: mpp_update_domains
303 
304  implicit none
305  type(fv3jedi_geom), intent(inout) :: geom
306  real(kind=kind_real), intent(inout) :: psi(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz) !Stream function
307  real(kind=kind_real), intent(inout) :: chi(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz) !Velocity potential
308  real(kind=kind_real), intent(out) :: ua(geom%isc:geom%iec,geom%jsc:geom%jec,1:geom%npz) !Agrid winds (u)
309  real(kind=kind_real), intent(out) :: va(geom%isc:geom%iec,geom%jsc:geom%jec,1:geom%npz) !Agrid winds (v)
310 
311  integer :: i,j,k
312 
313  ! x-----------------x
314  ! | |
315  ! | |
316  ! | |
317  ! | x |
318  ! | psi,chi |
319  ! | ua,va |
320  ! | |
321  ! | |
322  ! x-----------------x
323 
324  call mpp_update_domains(psi, geom%domain, complete=.true.)
325  call mpp_update_domains(chi, geom%domain, complete=.true.)
326 
327  do k=1,geom%npz
328  do j=geom%jsc,geom%jec
329  do i=geom%isc,geom%iec
330 
331  ua(i,j,k) = (psi(i,j+1,k) - psi(i,j-1,k))/(geom%dyc(i,j) + geom%dyc(i,j+1)) + &
332  (chi(i+1,j,k) - chi(i-1,j,k))/(geom%dxc(i,j) + geom%dxc(i+1,j))
333  va(i,j,k) = -(psi(i+1,j,k) - psi(i-1,j,k))/(geom%dxc(i,j) + geom%dxc(i+1,j)) + &
334  (chi(i,j+1,k) - chi(i,j-1,k))/(geom%dyc(i,j) + geom%dyc(i,j+1))
335 
336  enddo
337  enddo
338  enddo
339 
340 end subroutine psichi_to_uava
341 
342 !----------------------------------------------------------------------------
343 
344 subroutine psichi_to_uava_adm(geom,psi_ad,chi_ad,ua_ad,va_ad)
345 
346  use fv_mp_adm_mod, only: mpp_update_domains_adm
347 
348  implicit none
349  type(fv3jedi_geom), intent(inout) :: geom
350  real(kind=kind_real), intent(inout) :: psi_ad(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz) !Stream function
351  real(kind=kind_real), intent(inout) :: chi_ad(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz) !Velocity potential
352  real(kind=kind_real), intent(inout) :: ua_ad(geom%isc:geom%iec,geom%jsc:geom%jec,1:geom%npz) !Agrid winds (u)
353  real(kind=kind_real), intent(inout) :: va_ad(geom%isc:geom%iec,geom%jsc:geom%jec,1:geom%npz) !Agrid winds (v)
354 
355  integer :: i,j,k
356  real(kind=kind_real) :: temp1
357  real(kind=kind_real) :: temp2
358  real(kind=kind_real) :: temp3
359  real(kind=kind_real) :: temp4
360 
361  ! x-----------------x
362  ! | |
363  ! | |
364  ! | |
365  ! | x |
366  ! | psi,chi |
367  ! | ua,va |
368  ! | |
369  ! | |
370  ! x-----------------x
371 
372  do k=geom%npz,1,-1
373  do j=geom%jec,geom%jsc,-1
374  do i=geom%iec,geom%isc,-1
375 
376  temp1 = va_ad(i,j,k)/(geom%dyc(i,j)+geom%dyc(i,j+1))
377  temp2 = -(va_ad(i,j,k)/(geom%dxc(i,j)+geom%dxc(i+1,j)))
378  chi_ad(i,j+1,k) = chi_ad(i,j+1,k) + temp1
379  chi_ad(i,j-1,k) = chi_ad(i,j-1,k) - temp1
380  psi_ad(i+1,j,k) = psi_ad(i+1,j,k) + temp2
381  psi_ad(i-1,j,k) = psi_ad(i-1,j,k) - temp2
382  va_ad(i,j,k) = 0.0_8
383 
384  temp3 = ua_ad(i,j,k)/(geom%dyc(i,j)+geom%dyc(i,j+1))
385  temp4 = ua_ad(i,j,k)/(geom%dxc(i,j)+geom%dxc(i+1,j))
386  psi_ad(i,j+1,k) = psi_ad(i,j+1,k) + temp3
387  psi_ad(i,j-1,k) = psi_ad(i,j-1,k) - temp3
388  chi_ad(i+1,j,k) = chi_ad(i+1,j,k) + temp4
389  chi_ad(i-1,j,k) = chi_ad(i-1,j,k) - temp4
390  ua_ad(i,j,k) = 0.0_8
391 
392  end do
393  end do
394  end do
395 
396  call mpp_update_domains_adm(ua_ad, chi_ad, geom%domain, complete=.true.)
397  call mpp_update_domains_adm(va_ad, psi_ad, geom%domain, complete=.true.)
398 
399 end subroutine psichi_to_uava_adm
400 
401 !----------------------------------------------------------------------------
402 
403 subroutine psichi_to_udvd(geom,psi,chi,u,v)
404 
405  use mpp_domains_mod, only: mpp_update_domains
406 
407  implicit none
408  type(fv3jedi_geom), intent(inout) :: geom
409  real(kind=kind_real), intent(inout) :: psi(geom%isd:geom%ied ,geom%jsd:geom%jed ,1:geom%npz) !Stream function
410  real(kind=kind_real), intent(inout) :: chi(geom%isd:geom%ied ,geom%jsd:geom%jed ,1:geom%npz) !Velocity potential
411  real(kind=kind_real), intent(out) :: u(geom%isd:geom%ied ,geom%jsd:geom%jed+1,1:geom%npz) !Dgrid winds (u)
412  real(kind=kind_real), intent(out) :: v(geom%isd:geom%ied+1,geom%jsd:geom%jed ,1:geom%npz) !Dgrid winds (v)
413 
414  integer :: i,j,k
415  real(kind=kind_real) :: chib(geom%isd:geom%ied,geom%jsd:geom%jed,1:geom%npz)
416 
417  ! x-----------------x
418  ! | |
419  ! | |
420  ! | |
421  ! vd| x |
422  ! | psi,chi |
423  ! | |
424  ! | |
425  ! | |
426  ! x-----------------x
427  ! ud
428 
429  !Fill halos of psi and chi
430  call mpp_update_domains(psi, geom%domain, complete=.true.)
431  call mpp_update_domains(chi, geom%domain, complete=.true.)
432 
433  !Interpolate chi to the B grid
434  call a2b_ord4(chi, chib, geom, geom%npx, geom%npy, geom%isc, geom%iec, geom%jsc, geom%jec, geom%halo)
435 
436  do k=1,geom%npz
437  do j=geom%jsc,geom%jec
438  do i=geom%isc,geom%iec
439 
440  u(i,j,k) = (psi(i,j+1,k) - psi(i,j,k))/(geom%dyc(i,j)) + &
441  (chib(i+1,j,k) - chib(i,j,k))/(geom%dx (i,j))
442  v(i,j,k) = -(psi(i+1,j,k) - psi(i,j,k))/(geom%dxc(i,j)) + &
443  (chib(i,j+1,k) - chib(i,j,k))/(geom%dy (i,j))
444 
445  enddo
446  enddo
447  enddo
448 
449 endsubroutine psichi_to_udvd
450 
451 !----------------------------------------------------------------------------
452 
453 
454 !----------------------------------------------------------------------------
455 
456 subroutine d_to_ac(geom, u, v, ua, va, uc, vc)
457 
458  !Convert D-grid winds to A-grid and C-grid
459  !Also fills the halo of u and v on the D-grid
460 
461  use mpp_domains_mod, only: mpp_update_domains, mpp_get_boundary, dgrid_ne
462 
463  implicit none
464 
465  type(fv3jedi_geom), intent(inout) :: geom
466  real(kind=kind_real), intent(inout) :: u(geom%isd:geom%ied ,geom%jsd:geom%jed+1,1:geom%npz)
467  real(kind=kind_real), intent(inout) :: v(geom%isd:geom%ied+1,geom%jsd:geom%jed ,1:geom%npz)
468  real(kind=kind_real), intent(inout) :: ua(geom%isd:geom%ied ,geom%jsd:geom%jed ,1:geom%npz)
469  real(kind=kind_real), intent(inout) :: va(geom%isd:geom%ied ,geom%jsd:geom%jed ,1:geom%npz)
470  real(kind=kind_real), optional, intent(inout) :: uc(geom%isd:geom%ied+1,geom%jsd:geom%jed ,1:geom%npz)
471  real(kind=kind_real), optional, intent(inout) :: vc(geom%isd:geom%ied ,geom%jsd:geom%jed+1,1:geom%npz)
472 
473  real(kind=kind_real) :: wbuffer(geom%jsc:geom%jec,geom%npz)
474  real(kind=kind_real) :: sbuffer(geom%isc:geom%iec,geom%npz)
475  real(kind=kind_real) :: ebuffer(geom%jsc:geom%jec,geom%npz)
476  real(kind=kind_real) :: nbuffer(geom%isc:geom%iec,geom%npz)
477 
478  integer isc,iec,jsc,jec
479  integer isd,ied,jsd,jed
480  integer npz
481  integer i,j,k
482 
483  real(kind=kind_real) :: ut(geom%isd:geom%ied, geom%jsd:geom%jed)
484  real(kind=kind_real) :: vt(geom%isd:geom%ied, geom%jsd:geom%jed)
485 
486  real(kind=kind_real) :: uctemp(geom%isd:geom%ied+1,geom%jsd:geom%jed ,1:geom%npz)
487  real(kind=kind_real) :: vctemp(geom%isd:geom%ied ,geom%jsd:geom%jed+1,1:geom%npz)
488 
489  isc=geom%isc
490  iec=geom%iec
491  jsc=geom%jsc
492  jec=geom%jec
493  isd=geom%isd
494  ied=geom%ied
495  jsd=geom%jsd
496  jed=geom%jed
497  npz = geom%npz
498 
499  uctemp = 0.0_kind_real
500  vctemp = 0.0_kind_real
501 
502  !Fill edge of patch winds
503  call mpp_get_boundary(u, v, geom%domain, &
504  wbuffery=wbuffer, ebuffery=ebuffer, &
505  sbufferx=sbuffer, nbufferx=nbuffer, &
506  gridtype=dgrid_ne, complete=.true. )
507  do k=1,npz
508  do i=isc,iec
509  u(i,jec+1,k) = nbuffer(i,k)
510  enddo
511  do j=jsc,jec
512  v(iec+1,j,k) = ebuffer(j,k)
513  enddo
514  enddo
515 
516  !Fill halo on winds
517  call mpp_update_domains(u, v, geom%domain, gridtype=dgrid_ne, complete=.true.)
518 
519  !Convert to C and A grids, winds perpendicular to grid cells
520  do k=1,npz
521  call d2a2c_vect(geom, u(:,:,k), v(:,:,k), &
522  ua(:,:,k), va(:,:,k), uctemp(:,:,k), vctemp(:,:,k), &
523  ut, vt, .true.)
524  enddo
525 
526  if (present(uc)) uc = uctemp
527  if (present(vc)) vc = vctemp
528 
529 end subroutine d_to_ac
530 
531 !----------------------------------------------------------------------------
532 
533 subroutine d2a2c_vect(geom, u, v, ua, va, uc, vc, ut, vt, dord4)
534 
535  implicit none
536  type(fv3jedi_geom), intent(IN), target :: geom
537  logical, intent(in):: dord4
538  real(kind=kind_real), intent(in) :: u(geom%isd:geom%ied ,geom%jsd:geom%jed+1)
539  real(kind=kind_real), intent(in) :: v(geom%isd:geom%ied+1,geom%jsd:geom%jed)
540  real(kind=kind_real), intent(out), dimension(geom%isd:geom%ied+1,geom%jsd:geom%jed ):: uc
541  real(kind=kind_real), intent(out), dimension(geom%isd:geom%ied ,geom%jsd:geom%jed+1):: vc
542  real(kind=kind_real), intent(out), dimension(geom%isd:geom%ied ,geom%jsd:geom%jed ):: ua, va
543  real(kind=kind_real), intent(out), dimension(geom%isd:geom%ied ,geom%jsd:geom%jed ):: ut, vt
544 
545  ! Local
546  real(kind=kind_real), dimension(geom%isd:geom%ied,geom%jsd:geom%jed):: utmp, vtmp
547  integer npt, i, j, ifirst, ilast, id
548  integer :: npx, npy
549  integer :: is, ie, js, je
550  integer :: isd, ied, jsd, jed
551 
552  real(kind=kind_real), pointer, dimension(:,:,:) :: sin_sg
553  real(kind=kind_real), pointer, dimension(:,:) :: cosa_u, cosa_v, cosa_s
554  real(kind=kind_real), pointer, dimension(:,:) :: rsin_u, rsin_v, rsin2
555  real(kind=kind_real), pointer, dimension(:,:) :: dxa,dya
556 
557  real(kind=kind_real), parameter :: a1 = 0.5625
558  real(kind=kind_real), parameter :: a2 = -0.0625
559  real(kind=kind_real), parameter :: c1 = -2./14.
560  real(kind=kind_real), parameter :: c2 = 11./14.
561  real(kind=kind_real), parameter :: c3 = 5./14.
562 
563  is = geom%isc
564  ie = geom%iec
565  js = geom%jsc
566  je = geom%jec
567  isd = geom%isd
568  ied = geom%ied
569  jsd = geom%jsd
570  jed = geom%jed
571  npx = geom%npx
572  npy = geom%npy
573 
574  sin_sg => geom%sin_sg
575  cosa_u => geom%cosa_u
576  cosa_v => geom%cosa_v
577  cosa_s => geom%cosa_s
578  rsin_u => geom%rsin_u
579  rsin_v => geom%rsin_v
580  rsin2 => geom%rsin2
581  dxa => geom%dxa
582  dya => geom%dya
583 
584  if ( dord4 ) then
585  id = 1
586  else
587  id = 0
588  endif
589 
590  npt = 4
591 
592 ! Initialize the non-existing corner regions
593  utmp(:,:) = 1.0e8_kind_real
594  vtmp(:,:) = 1.0e8_kind_real
595 
596  !----------
597  ! Interior:
598  !----------
599  do j=max(npt,js-1),min(npy-npt,je+1)
600  do i=max(npt,isd),min(npx-npt,ied)
601  utmp(i,j) = a2*(u(i,j-1)+u(i,j+2)) + a1*(u(i,j)+u(i,j+1))
602  enddo
603  enddo
604  do j=max(npt,jsd),min(npy-npt,jed)
605  do i=max(npt,is-1),min(npx-npt,ie+1)
606  vtmp(i,j) = a2*(v(i-1,j)+v(i+2,j)) + a1*(v(i,j)+v(i+1,j))
607  enddo
608  enddo
609 
610  !----------
611  ! edges:
612  !----------
613  if ( js==1 .or. jsd<npt) then
614  do j=jsd,npt-1
615  do i=isd,ied
616  utmp(i,j) = 0.5*(u(i,j) + u(i,j+1))
617  vtmp(i,j) = 0.5*(v(i,j) + v(i+1,j))
618  enddo
619  enddo
620  endif
621  if ( (je+1)==npy .or. jed>=(npy-npt)) then
622  do j=npy-npt+1,jed
623  do i=isd,ied
624  utmp(i,j) = 0.5*(u(i,j) + u(i,j+1))
625  vtmp(i,j) = 0.5*(v(i,j) + v(i+1,j))
626  enddo
627  enddo
628  endif
629 
630  if ( is==1 .or. isd<npt ) then
631  do j=max(npt,jsd),min(npy-npt,jed)
632  do i=isd,npt-1
633  utmp(i,j) = 0.5*(u(i,j) + u(i,j+1))
634  vtmp(i,j) = 0.5*(v(i,j) + v(i+1,j))
635  enddo
636  enddo
637  endif
638  if ( (ie+1)==npx .or. ied>=(npx-npt)) then
639  do j=max(npt,jsd),min(npy-npt,jed)
640  do i=npx-npt+1,ied
641  utmp(i,j) = 0.5*(u(i,j) + u(i,j+1))
642  vtmp(i,j) = 0.5*(v(i,j) + v(i+1,j))
643  enddo
644  enddo
645  endif
646 
647 ! Contra-variant components at cell center:
648  do j=js-1-id,je+1+id
649  do i=is-1-id,ie+1+id
650  ua(i,j) = (utmp(i,j)-vtmp(i,j)*cosa_s(i,j)) * rsin2(i,j)
651  va(i,j) = (vtmp(i,j)-utmp(i,j)*cosa_s(i,j)) * rsin2(i,j)
652  enddo
653  enddo
654 
655 ! A -> C
656 !--------------
657 ! Fix the edges
658 !--------------
659 ! Xdir:
660  if( geom%sw_corner ) then
661  do i=-2,0
662  utmp(i,0) = -vtmp(0,1-i)
663  enddo
664  endif
665  if( geom%se_corner ) then
666  do i=0,2
667  utmp(npx+i,0) = vtmp(npx,i+1)
668  enddo
669  endif
670  if( geom%ne_corner ) then
671  do i=0,2
672  utmp(npx+i,npy) = -vtmp(npx,je-i)
673  enddo
674  endif
675  if( geom%nw_corner ) then
676  do i=-2,0
677  utmp(i,npy) = vtmp(0,je+i)
678  enddo
679  endif
680 
681  ifirst = max(3, is-1)
682  ilast = min(npx-2,ie+2)
683 
684 !---------------------------------------------
685 ! 4th order interpolation for interior points:
686 !---------------------------------------------
687  do j=js-1,je+1
688  do i=ifirst,ilast
689  uc(i,j) = a2*(utmp(i-2,j)+utmp(i+1,j)) + a1*(utmp(i-1,j)+utmp(i,j))
690  ut(i,j) = (uc(i,j) - v(i,j)*cosa_u(i,j))*rsin_u(i,j)
691  enddo
692  enddo
693 
694 ! Xdir:
695  if( geom%sw_corner ) then
696  ua(-1,0) = -va(0,2)
697  ua( 0,0) = -va(0,1)
698  endif
699  if( geom%se_corner ) then
700  ua(npx, 0) = va(npx,1)
701  ua(npx+1,0) = va(npx,2)
702  endif
703  if( geom%ne_corner ) then
704  ua(npx, npy) = -va(npx,npy-1)
705  ua(npx+1,npy) = -va(npx,npy-2)
706  endif
707  if( geom%nw_corner ) then
708  ua(-1,npy) = va(0,npy-2)
709  ua( 0,npy) = va(0,npy-1)
710  endif
711 
712  if( is==1 ) then
713  do j=js-1,je+1
714  uc(0,j) = c1*utmp(-2,j) + c2*utmp(-1,j) + c3*utmp(0,j)
715  ut(1,j) = edge_interpolate4(ua(-1:2,j), dxa(-1:2,j))
716  !Want to use the UPSTREAM value
717  if (ut(1,j) > 0.) then
718  uc(1,j) = ut(1,j)*sin_sg(0,j,3)
719  else
720  uc(1,j) = ut(1,j)*sin_sg(1,j,1)
721  end if
722  uc(2,j) = c1*utmp(3,j) + c2*utmp(2,j) + c3*utmp(1,j)
723  ut(0,j) = (uc(0,j) - v(0,j)*cosa_u(0,j))*rsin_u(0,j)
724  ut(2,j) = (uc(2,j) - v(2,j)*cosa_u(2,j))*rsin_u(2,j)
725  enddo
726  endif
727 
728  if( (ie+1)==npx ) then
729  do j=js-1,je+1
730  uc(npx-1,j) = c1*utmp(npx-3,j)+c2*utmp(npx-2,j)+c3*utmp(npx-1,j)
731  ut(npx, j) = edge_interpolate4(ua(npx-2:npx+1,j), dxa(npx-2:npx+1,j))
732  if (ut(npx,j) > 0.) then
733  uc(npx,j) = ut(npx,j)*sin_sg(npx-1,j,3)
734  else
735  uc(npx,j) = ut(npx,j)*sin_sg(npx,j,1)
736  end if
737  uc(npx+1,j) = c3*utmp(npx,j) + c2*utmp(npx+1,j) + c1*utmp(npx+2,j)
738  ut(npx-1,j) = (uc(npx-1,j)-v(npx-1,j)*cosa_u(npx-1,j))*rsin_u(npx-1,j)
739  ut(npx+1,j) = (uc(npx+1,j)-v(npx+1,j)*cosa_u(npx+1,j))*rsin_u(npx+1,j)
740  enddo
741  endif
742 
743 
744 !------
745 ! Ydir:
746 !------
747  if( geom%sw_corner ) then
748  do j=-2,0
749  vtmp(0,j) = -utmp(1-j,0)
750  enddo
751  endif
752  if( geom%nw_corner ) then
753  do j=0,2
754  vtmp(0,npy+j) = utmp(j+1,npy)
755  enddo
756  endif
757  if( geom%se_corner ) then
758  do j=-2,0
759  vtmp(npx,j) = utmp(ie+j,0)
760  enddo
761  endif
762  if( geom%ne_corner ) then
763  do j=0,2
764  vtmp(npx,npy+j) = -utmp(ie-j,npy)
765  enddo
766  endif
767  if( geom%sw_corner ) then
768  va(0,-1) = -ua(2,0)
769  va(0, 0) = -ua(1,0)
770  endif
771  if( geom%se_corner ) then
772  va(npx, 0) = ua(npx-1,0)
773  va(npx,-1) = ua(npx-2,0)
774  endif
775  if( geom%ne_corner ) then
776  va(npx,npy ) = -ua(npx-1,npy)
777  va(npx,npy+1) = -ua(npx-2,npy)
778  endif
779  if( geom%nw_corner ) then
780  va(0,npy) = ua(1,npy)
781  va(0,npy+1) = ua(2,npy)
782  endif
783 
784 
785  do j=js-1,je+2
786  if ( j==1) then
787  do i=is-1,ie+1
788  vt(i,j) = edge_interpolate4(va(i,-1:2), dya(i,-1:2))
789  if (vt(i,j) > 0.) then
790  vc(i,j) = vt(i,j)*sin_sg(i,j-1,4)
791  else
792  vc(i,j) = vt(i,j)*sin_sg(i,j,2)
793  end if
794  enddo
795  elseif ( j==0 .or. j==(npy-1)) then
796  do i=is-1,ie+1
797  vc(i,j) = c1*vtmp(i,j-2) + c2*vtmp(i,j-1) + c3*vtmp(i,j)
798  vt(i,j) = (vc(i,j) - u(i,j)*cosa_v(i,j))*rsin_v(i,j)
799  enddo
800  elseif ( j==2 .or. j==(npy+1)) then
801  do i=is-1,ie+1
802  vc(i,j) = c1*vtmp(i,j+1) + c2*vtmp(i,j) + c3*vtmp(i,j-1)
803  vt(i,j) = (vc(i,j) - u(i,j)*cosa_v(i,j))*rsin_v(i,j)
804  enddo
805  elseif ( j==npy ) then
806  do i=is-1,ie+1
807  vt(i,j) = edge_interpolate4(va(i,j-2:j+1), dya(i,j-2:j+1))
808  if (vt(i,j) > 0.) then
809  vc(i,j) = vt(i,j)*sin_sg(i,j-1,4)
810  else
811  vc(i,j) = vt(i,j)*sin_sg(i,j,2)
812  end if
813  enddo
814  else
815 ! 4th order interpolation for interior points:
816  do i=is-1,ie+1
817  vc(i,j) = a2*(vtmp(i,j-2)+vtmp(i,j+1)) + a1*(vtmp(i,j-1)+vtmp(i,j))
818  vt(i,j) = (vc(i,j) - u(i,j)*cosa_v(i,j))*rsin_v(i,j)
819  enddo
820  endif
821  enddo
822 
823  end subroutine d2a2c_vect
824 
825 !----------------------------------------------------------------------------
826 
827 real(kind=kind_real) function edge_interpolate4(ua, dxa)
828 
829  implicit none
830 
831  real(kind=kind_real), intent(in) :: ua(4)
832  real(kind=kind_real), intent(in) :: dxa(4)
833 
834  real(kind=kind_real) :: t1, t2
835 
836  t1 = dxa(1) + dxa(2)
837  t2 = dxa(3) + dxa(4)
838  edge_interpolate4 = 0.5*( ((t1+dxa(2))*ua(2)-dxa(2)*ua(1)) / t1 + &
839  ((t2+dxa(3))*ua(3)-dxa(3)*ua(4)) / t2 )
840 
841 end function edge_interpolate4
842 
843 !----------------------------------------------------------------------------
844 
845  subroutine divergence_corner(geom, u, v, ua, va, divg_d)
846 
847  !> Refactored from sw_core.F90
848 
849  !> Divergence on the B-Grid from D and A grid winds
850 
851  implicit none
852  type(fv3jedi_geom), intent(in) :: geom
853  real(kind=kind_real), intent(in), dimension(geom%isd:geom%ied, geom%jsd:geom%jed+1):: u
854  real(kind=kind_real), intent(in), dimension(geom%isd:geom%ied+1,geom%jsd:geom%jed ):: v
855  real(kind=kind_real), intent(in), dimension(geom%isd:geom%ied ,geom%jsd:geom%jed ):: ua
856  real(kind=kind_real), intent(in), dimension(geom%isd:geom%ied ,geom%jsd:geom%jed ):: va
857  real(kind=kind_real), intent(out), dimension(geom%isd:geom%ied+1,geom%jsd:geom%jed+1):: divg_d
858 
859  !Local
860  integer :: is, ie, js, je
861  integer :: npx, npy
862  integer i,j,is2,ie1
863  real(kind=kind_real) :: uf(geom%isc-2:geom%iec+2,geom%jsc-1:geom%jec+2)
864  real(kind=kind_real) :: vf(geom%isc-1:geom%iec+2,geom%jsc-2:geom%jec+2)
865 
866  is = geom%isc
867  ie = geom%iec
868  js = geom%jsc
869  je = geom%jec
870 
871  npx = geom%npx
872  npy = geom%npy
873 
874  is2 = max(2,is)
875  ie1 = min(npx-1,ie+1)
876 
877  !uf
878  do j=js,je+1
879  if ( j==1 .or. j==npy ) then
880  do i=is-1,ie+1
881  uf(i,j) = u(i,j)*geom%dyc(i,j)*0.5*(geom%sin_sg(i,j-1,4)+geom%sin_sg(i,j,2))
882  enddo
883  else
884  do i=is-1,ie+1
885  uf(i,j) = (u(i,j)-0.25*(va(i,j-1)+va(i,j))*(geom%cos_sg(i,j-1,4)+geom%cos_sg(i,j,2))) &
886  * geom%dyc(i,j)*0.5*(geom%sin_sg(i,j-1,4)+geom%sin_sg(i,j,2))
887  enddo
888  endif
889  enddo
890 
891  !vf
892  do j=js-1,je+1
893  do i=is2,ie1
894  vf(i,j) = (v(i,j) - 0.25*(ua(i-1,j)+ua(i,j))*(geom%cos_sg(i-1,j,3)+geom%cos_sg(i,j,1))) &
895  *geom%dxc(i,j)*0.5*(geom%sin_sg(i-1,j,3)+geom%sin_sg(i,j,1))
896  enddo
897  if ( is == 1 ) vf(1, j) = v(1, j)*geom%dxc(1, j)*0.5*(geom%sin_sg(0,j,3)+geom%sin_sg(1,j,1))
898  if ( (ie+1)==npx ) vf(npx,j) = v(npx,j)*geom%dxc(npx,j)*0.5*(geom%sin_sg(npx-1,j,3)+geom%sin_sg(npx,j,1))
899  enddo
900 
901  !divg_d 1
902  do j=js,je+1
903  do i=is,ie+1
904  divg_d(i,j) = (vf(i,j-1) - vf(i,j)) + (uf(i-1,j) - uf(i,j))
905  enddo
906  enddo
907 
908  !Remove the extra term at the corners:
909  if (geom%sw_corner) divg_d(1, 1) = divg_d(1, 1) - vf(1, 0)
910  if (geom%se_corner) divg_d(npx, 1) = divg_d(npx, 1) - vf(npx, 0)
911  if (geom%ne_corner) divg_d(npx,npy) = divg_d(npx,npy) + vf(npx,npy)
912  if (geom%nw_corner) divg_d(1, npy) = divg_d(1, npy) + vf(1, npy)
913 
914  !Mult by grid area
915  do j=js,je+1
916  do i=is,ie+1
917  divg_d(i,j) = geom%rarea_c(i,j)*divg_d(i,j)
918  enddo
919  enddo
920 
921 end subroutine divergence_corner
922 
923 !----------------------------------------------------------------------------
924 
925  subroutine a2b_ord4(qin, qout, geom, npx, npy, is, ie, js, je, ng)
926 
927  implicit none
928 
929  integer, intent(IN):: npx, npy, is, ie, js, je, ng
930  real(kind=kind_real), intent(IN):: qin(is-ng:ie+ng,js-ng:je+ng) ! A-grid field
931  real(kind=kind_real), intent(OUT):: qout(is-ng:ie+ng,js-ng:je+ng) ! Output B-grid field
932  type(fv3jedi_geom), intent(IN), target :: geom
933 
934  real(kind=kind_real) qx(is:ie+1,js-ng:je+ng)
935  real(kind=kind_real) qy(is-ng:ie+ng,js:je+1)
936  real(kind=kind_real) qxx(is-ng:ie+ng,js-ng:je+ng)
937  real(kind=kind_real) qyy(is-ng:ie+ng,js-ng:je+ng)
938  real(kind=kind_real) g_in, g_ou
939  real(kind=kind_real):: p0(2)
940  real(kind=kind_real):: q1(is-1:ie+1), q2(js-1:je+1)
941  integer:: i, j, is1, js1, is2, js2, ie1, je1
942 
943  real(kind=kind_real), parameter :: c1 = 2./3.
944  real(kind=kind_real), parameter :: c2 = -1./6.
945  real(kind=kind_real), parameter :: r3 = 1./3.
946  real(kind=kind_real), parameter :: a1 = 0.5625 ! 9/16
947  real(kind=kind_real), parameter :: a2 = -0.0625 ! -1/16
948  real(kind=kind_real), parameter :: b1 = 7./12. ! 0.58333333
949  real(kind=kind_real), parameter :: b2 = -1./12.
950 
951  is1 = max(1,is-1)
952  js1 = max(1,js-1)
953  is2 = max(2,is)
954  js2 = max(2,js)
955 
956  ie1 = min(npx-1,ie+1)
957  je1 = min(npy-1,je+1)
958 
959 ! Corners:
960 ! 3-way extrapolation
961 
962  if ( geom%sw_corner ) then
963  p0(1:2) = geom%grid(1,1,1:2)
964  qout(1,1) = (extrap_corner(p0, geom%agrid(1,1,1:2), geom%agrid( 2, 2,1:2), qin(1,1), qin( 2, 2)) + &
965  extrap_corner(p0, geom%agrid(0,1,1:2), geom%agrid(-1, 2,1:2), qin(0,1), qin(-1, 2)) + &
966  extrap_corner(p0, geom%agrid(1,0,1:2), geom%agrid( 2,-1,1:2), qin(1,0), qin( 2,-1)))*r3
967 
968  endif
969  if ( geom%se_corner ) then
970  p0(1:2) = geom%grid(npx,1,1:2)
971  qout(npx,1) = (extrap_corner(p0, geom%agrid(npx-1,1,1:2), geom%agrid(npx-2, 2,1:2), qin(npx-1,1), qin(npx-2, 2)) + &
972  extrap_corner(p0, geom%agrid(npx-1,0,1:2), geom%agrid(npx-2,-1,1:2), qin(npx-1,0), qin(npx-2,-1)) + &
973  extrap_corner(p0, geom%agrid(npx ,1,1:2), geom%agrid(npx+1, 2,1:2), qin(npx ,1), qin(npx+1, 2)))*r3
974  endif
975  if ( geom%ne_corner ) then
976  p0(1:2) = geom%grid(npx,npy,1:2)
977  qout(npx,npy) = (extrap_corner(p0, geom%agrid(npx-1,npy-1,1:2), geom%agrid(npx-2,npy-2,1:2), qin(npx-1,npy-1), qin(npx-2,npy-2)) + &
978  extrap_corner(p0, geom%agrid(npx ,npy-1,1:2), geom%agrid(npx+1,npy-2,1:2), qin(npx ,npy-1), qin(npx+1,npy-2)) + &
979  extrap_corner(p0, geom%agrid(npx-1,npy ,1:2), geom%agrid(npx-2,npy+1,1:2), qin(npx-1,npy ), qin(npx-2,npy+1)))*r3
980  endif
981  if ( geom%nw_corner ) then
982  p0(1:2) = geom%grid(1,npy,1:2)
983  qout(1,npy) = (extrap_corner(p0, geom%agrid(1,npy-1,1:2), geom%agrid( 2,npy-2,1:2), qin(1,npy-1), qin( 2,npy-2)) + &
984  extrap_corner(p0, geom%agrid(0,npy-1,1:2), geom%agrid(-1,npy-2,1:2), qin(0,npy-1), qin(-1,npy-2)) + &
985  extrap_corner(p0, geom%agrid(1,npy, 1:2), geom%agrid( 2,npy+1,1:2), qin(1,npy ), qin( 2,npy+1)))*r3
986  endif
987 
988 !------------
989 ! X-Interior:
990 !------------
991  do j=max(1,js-2),min(npy-1,je+2)
992  do i=max(3,is), min(npx-2,ie+1)
993  qx(i,j) = b2*(qin(i-2,j)+qin(i+1,j)) + b1*(qin(i-1,j)+qin(i,j))
994  enddo
995  enddo
996 
997  ! *** West Edges:
998  if ( is==1 ) then
999  do j=js1, je1
1000  q2(j) = (qin(0,j)*geom%dxa(1,j) + qin(1,j)*geom%dxa(0,j))/(geom%dxa(0,j) + geom%dxa(1,j))
1001  enddo
1002  do j=js2, je1
1003  qout(1,j) = geom%edge_w(j)*q2(j-1) + (1.-geom%edge_w(j))*q2(j)
1004  enddo
1005 !
1006  do j=max(1,js-2),min(npy-1,je+2)
1007  g_in = geom%dxa(2,j) / geom%dxa(1,j)
1008  g_ou = geom%dxa(-1,j) / geom%dxa(0,j)
1009  qx(1,j) = 0.5*( ((2.+g_in)*qin(1,j)-qin( 2,j))/(1.+g_in) + &
1010  ((2.+g_ou)*qin(0,j)-qin(-1,j))/(1.+g_ou) )
1011  qx(2,j) = ( 3.*(g_in*qin(1,j)+qin(2,j))-(g_in*qx(1,j)+qx(3,j)) ) / (2.+2.*g_in)
1012  enddo
1013  endif
1014 
1015  ! East Edges:
1016  if ( (ie+1)==npx ) then
1017  do j=js1, je1
1018  q2(j) = (qin(npx-1,j)*geom%dxa(npx,j) + qin(npx,j)*geom%dxa(npx-1,j))/(geom%dxa(npx-1,j) + geom%dxa(npx,j))
1019  enddo
1020  do j=js2, je1
1021  qout(npx,j) = geom%edge_e(j)*q2(j-1) + (1.-geom%edge_e(j))*q2(j)
1022  enddo
1023 !
1024  do j=max(1,js-2),min(npy-1,je+2)
1025  g_in = geom%dxa(npx-2,j) / geom%dxa(npx-1,j)
1026  g_ou = geom%dxa(npx+1,j) / geom%dxa(npx,j)
1027  qx(npx,j) = 0.5*( ((2.+g_in)*qin(npx-1,j)-qin(npx-2,j))/(1.+g_in) + &
1028  ((2.+g_ou)*qin(npx, j)-qin(npx+1,j))/(1.+g_ou) )
1029  qx(npx-1,j) = (3.*(qin(npx-2,j)+g_in*qin(npx-1,j)) - (g_in*qx(npx,j)+qx(npx-2,j)))/(2.+2.*g_in)
1030  enddo
1031  endif
1032 
1033 !------------
1034 ! Y-Interior:
1035 !------------
1036 
1037  do j=max(3,js),min(npy-2,je+1)
1038  do i=max(1,is-2), min(npx-1,ie+2)
1039  qy(i,j) = b2*(qin(i,j-2)+qin(i,j+1)) + b1*(qin(i,j-1) + qin(i,j))
1040  enddo
1041  enddo
1042 
1043  ! South Edges:
1044  if ( js==1 ) then
1045  do i=is1, ie1
1046  q1(i) = (qin(i,0)*geom%dya(i,1) + qin(i,1)*geom%dya(i,0))/(geom%dya(i,0) + geom%dya(i,1))
1047  enddo
1048  do i=is2, ie1
1049  qout(i,1) = geom%edge_s(i)*q1(i-1) + (1.-geom%edge_s(i))*q1(i)
1050  enddo
1051 !
1052  do i=max(1,is-2),min(npx-1,ie+2)
1053  g_in = geom%dya(i,2) / geom%dya(i,1)
1054  g_ou = geom%dya(i,-1) / geom%dya(i,0)
1055  qy(i,1) = 0.5*( ((2.+g_in)*qin(i,1)-qin(i,2))/(1.+g_in) + &
1056  ((2.+g_ou)*qin(i,0)-qin(i,-1))/(1.+g_ou) )
1057  qy(i,2) = (3.*(g_in*qin(i,1)+qin(i,2)) - (g_in*qy(i,1)+qy(i,3)))/(2.+2.*g_in)
1058  enddo
1059  endif
1060 
1061  ! North Edges:
1062  if ( (je+1)==npy ) then
1063  do i=is1, ie1
1064  q1(i) = (qin(i,npy-1)*geom%dya(i,npy) + qin(i,npy)*geom%dya(i,npy-1))/(geom%dya(i,npy-1)+geom%dya(i,npy))
1065  enddo
1066  do i=is2, ie1
1067  qout(i,npy) = geom%edge_n(i)*q1(i-1) + (1.-geom%edge_n(i))*q1(i)
1068  enddo
1069 !
1070  do i=max(1,is-2),min(npx-1,ie+2)
1071  g_in = geom%dya(i,npy-2) / geom%dya(i,npy-1)
1072  g_ou = geom%dya(i,npy+1) / geom%dya(i,npy)
1073  qy(i,npy) = 0.5*( ((2.+g_in)*qin(i,npy-1)-qin(i,npy-2))/(1.+g_in) + &
1074  ((2.+g_ou)*qin(i,npy )-qin(i,npy+1))/(1.+g_ou) )
1075  qy(i,npy-1) = (3.*(qin(i,npy-2)+g_in*qin(i,npy-1)) - (g_in*qy(i,npy)+qy(i,npy-2)))/(2.+2.*g_in)
1076  enddo
1077  endif
1078 
1079 !--------------------------------------
1080 
1081  do j=max(3,js),min(npy-2,je+1)
1082  do i=max(2,is),min(npx-1,ie+1)
1083  qxx(i,j) = a2*(qx(i,j-2)+qx(i,j+1)) + a1*(qx(i,j-1)+qx(i,j))
1084  enddo
1085  enddo
1086 
1087  if ( js==1 ) then
1088  do i=max(2,is),min(npx-1,ie+1)
1089  qxx(i,2) = c1*(qx(i,1)+qx(i,2))+c2*(qout(i,1)+qxx(i,3))
1090  enddo
1091  endif
1092  if ( (je+1)==npy ) then
1093  do i=max(2,is),min(npx-1,ie+1)
1094  qxx(i,npy-1) = c1*(qx(i,npy-2)+qx(i,npy-1))+c2*(qout(i,npy)+qxx(i,npy-2))
1095  enddo
1096  endif
1097 
1098 
1099  do j=max(2,js),min(npy-1,je+1)
1100  do i=max(3,is),min(npx-2,ie+1)
1101  qyy(i,j) = a2*(qy(i-2,j)+qy(i+1,j)) + a1*(qy(i-1,j)+qy(i,j))
1102  enddo
1103  if ( is==1 ) qyy(2,j) = c1*(qy(1,j)+qy(2,j))+c2*(qout(1,j)+qyy(3,j))
1104  if((ie+1)==npx) qyy(npx-1,j) = c1*(qy(npx-2,j)+qy(npx-1,j))+c2*(qout(npx,j)+qyy(npx-2,j))
1105 
1106  do i=max(2,is),min(npx-1,ie+1)
1107  qout(i,j) = 0.5*(qxx(i,j) + qyy(i,j)) ! averaging
1108  enddo
1109  enddo
1110 
1111  end subroutine a2b_ord4
1112 
1113 !----------------------------------------------------------------------------
1114 
1115 real(kind=kind_real) function extrap_corner ( p0, p1, p2, q1, q2 )
1116 
1117  implicit none
1118  real(kind=kind_real), intent(in ), dimension(2) :: p0, p1, p2
1119  real(kind=kind_real), intent(in ) :: q1, q2
1120  real(kind=kind_real) :: x1, x2
1121 
1122  x1 = great_circle_dist( p1, p0 )
1123  x2 = great_circle_dist( p2, p0 )
1124 
1125  extrap_corner = q1 + x1/(x2-x1) * (q1-q2)
1126 
1127  end function extrap_corner
1128 
1129 !----------------------------------------------------------------------------
1130 
1131  real(kind=kind_real) function great_circle_dist( q1, q2 )
1132 
1133  implicit none
1134  real(kind=kind_real), intent(IN) :: q1(2), q2(2)
1135 
1136  real (kind=kind_real):: p1(2), p2(2)
1137  integer n
1138 
1139  do n=1,2
1140  p1(n) = q1(n)
1141  p2(n) = q2(n)
1142  enddo
1143 
1144  great_circle_dist = asin( sqrt( sin((p1(2)-p2(2))/2.)**2 + cos(p1(2))*cos(p2(2))* &
1145  sin((p1(1)-p2(1))/2.)**2 ) ) * 2.
1146 
1147  end function great_circle_dist
1148 
1149 !----------------------------------------------------------------------------
1150 
1151  SUBROUTINE psichi_to_udvd_adm(geom, psi_ad, chi_ad, u_ad, v_ad)
1152 
1153  use fv_mp_adm_mod, only: mpp_update_domains_adm
1154 
1155  IMPLICIT NONE
1156  TYPE(fv3jedi_geom), INTENT(INOUT) :: geom
1157 !Stream function
1158  REAL(kind=kind_real) :: psi(geom%isd:geom%ied, &
1159 & geom%jsd:geom%jed, geom%npz)
1160  REAL(kind=kind_real), INTENT(INOUT) :: psi_ad(geom%isd:geom%&
1161 & ied, geom%jsd:geom%jed, geom%npz)
1162 !Velocity potential
1163  REAL(kind=kind_real) :: chi(geom%isd:geom%ied, &
1164 & geom%jsd:geom%jed, geom%npz)
1165  REAL(kind=kind_real), INTENT(INOUT) :: chi_ad(geom%isd:geom%&
1166 & ied, geom%jsd:geom%jed, geom%npz)
1167 !Dgrid winds (u)
1168  REAL(kind=kind_real) :: u(geom%isd:geom%ied, geom%jsd:geom%&
1169 & jed+1, geom%npz)
1170  REAL(kind=kind_real) :: u_ad(geom%isd:geom%ied, geom%jsd:&
1171 & geom%jed+1, geom%npz)
1172 !Dgrid winds (v)
1173  REAL(kind=kind_real) :: v(geom%isd:geom%ied+1, geom%jsd:&
1174 & geom%jed, geom%npz)
1175  REAL(kind=kind_real) :: v_ad(geom%isd:geom%ied+1, geom%jsd:&
1176 & geom%jed, geom%npz)
1177  INTEGER :: i, j, k
1178  REAL(kind=kind_real) :: chib(geom%isd:geom%ied, geom%jsd:&
1179 & geom%jed, geom%npz)
1180  REAL(kind=kind_real) :: chib_ad(geom%isd:geom%ied, geom%jsd&
1181 & :geom%jed, geom%npz)
1182 ! x-----------------x
1183 ! | |
1184 ! | |
1185 ! | |
1186 ! vd| x |
1187 ! | psi,chi |
1188 ! | |
1189 ! | |
1190 ! | |
1191 ! x-----------------x
1192 ! ud
1193 !Fill halos of psi and chi
1194 !Interpolate chi to the B grid
1195  REAL(kind=kind_real) :: temp_ad
1196  REAL(kind=kind_real) :: temp_ad0
1197  REAL(kind=kind_real) :: temp_ad1
1198  REAL(kind=kind_real) :: temp_ad2
1199 
1200  u = 0.0_kind_real
1201  v = 0.0_kind_real
1202  psi = 0.0_kind_real
1203  chi = 0.0_kind_real
1204 
1205  chib_ad = 0.0_8
1206  DO k=geom%npz,1,-1
1207  DO j=geom%jec,geom%jsc,-1
1208  DO i=geom%iec,geom%isc,-1
1209  temp_ad = v_ad(i, j, k)/geom%dy(i, j)
1210  temp_ad0 = -(v_ad(i, j, k)/geom%dxc(i, j))
1211  chib_ad(i, j+1, k) = chib_ad(i, j+1, k) + temp_ad
1212  chib_ad(i, j, k) = chib_ad(i, j, k) - temp_ad
1213  psi_ad(i+1, j, k) = psi_ad(i+1, j, k) + temp_ad0
1214  psi_ad(i, j, k) = psi_ad(i, j, k) - temp_ad0
1215  v_ad(i, j, k) = 0.0_8
1216  temp_ad1 = u_ad(i, j, k)/geom%dyc(i, j)
1217  temp_ad2 = u_ad(i, j, k)/geom%dx(i, j)
1218  psi_ad(i, j+1, k) = psi_ad(i, j+1, k) + temp_ad1
1219  psi_ad(i, j, k) = psi_ad(i, j, k) - temp_ad1
1220  chib_ad(i+1, j, k) = chib_ad(i+1, j, k) + temp_ad2
1221  chib_ad(i, j, k) = chib_ad(i, j, k) - temp_ad2
1222  u_ad(i, j, k) = 0.0_8
1223  END DO
1224  END DO
1225  END DO
1226  CALL a2b_ord4_adm(chi, chi_ad, chib, chib_ad, geom, geom%npx, geom%&
1227 & npy, geom%isc, geom%iec, geom%jsc, geom%jec&
1228 & , geom%halo)
1229  CALL mpp_update_domains_adm(chi, chi_ad, geom%domain, complete=&
1230 & .true.)
1231  CALL mpp_update_domains_adm(psi, psi_ad, geom%domain, complete=&
1232 & .true.)
1233  END SUBROUTINE psichi_to_udvd_adm
1234 
1235 !----------------------------------------------------------------------------
1236 
1237  SUBROUTINE a2b_ord4_adm(qin, qin_ad, qout, qout_ad, geom, npx, npy, is&
1238 & , ie, js, je, ng)
1239  IMPLICIT NONE
1240  INTEGER, INTENT(IN) :: npx, npy, is, ie, js, je, ng
1241 ! A-grid field
1242  REAL(kind=kind_real), INTENT(IN) :: qin(is-ng:ie+ng, js-ng:je+ng)
1243  REAL(kind=kind_real) :: qin_ad(is-ng:ie+ng, js-ng:je+ng)
1244 ! Output B-grid field
1245  REAL(kind=kind_real) :: qout(is-ng:ie+ng, js-ng:je+ng)
1246  REAL(kind=kind_real) :: qout_ad(is-ng:ie+ng, js-ng:je+ng)
1247  TYPE(fv3jedi_geom), INTENT(IN), TARGET :: geom
1248  REAL(kind=kind_real) :: qx(is:ie+1, js-ng:je+ng)
1249  REAL(kind=kind_real) :: qx_ad(is:ie+1, js-ng:je+ng)
1250  REAL(kind=kind_real) :: qy(is-ng:ie+ng, js:je+1)
1251  REAL(kind=kind_real) :: qy_ad(is-ng:ie+ng, js:je+1)
1252  REAL(kind=kind_real) :: qxx(is-ng:ie+ng, js-ng:je+ng)
1253  REAL(kind=kind_real) :: qxx_ad(is-ng:ie+ng, js-ng:je+ng)
1254  REAL(kind=kind_real) :: qyy(is-ng:ie+ng, js-ng:je+ng)
1255  REAL(kind=kind_real) :: qyy_ad(is-ng:ie+ng, js-ng:je+ng)
1256  REAL(kind=kind_real) :: g_in, g_ou
1257  REAL(kind=kind_real) :: p0(2)
1258  REAL(kind=kind_real) :: q1(is-1:ie+1), q2(js-1:je+1)
1259  REAL(kind=kind_real) :: q1_ad(is-1:ie+1), q2_ad(js-1:je+1)
1260  INTEGER :: i, j, is1, js1, is2, js2, ie1, je1
1261  REAL(kind=kind_real), PARAMETER :: c1=2./3.
1262  REAL(kind=kind_real), PARAMETER :: c2=-(1./6.)
1263  REAL(kind=kind_real), PARAMETER :: r3=1./3.
1264 ! 9/16
1265  REAL(kind=kind_real), PARAMETER :: a1=0.5625
1266 ! -1/16
1267  REAL(kind=kind_real), PARAMETER :: a2=-0.0625
1268 ! 0.58333333
1269  REAL(kind=kind_real), PARAMETER :: b1=7./12.
1270  REAL(kind=kind_real), PARAMETER :: b2=-(1./12.)
1271  INTRINSIC max
1272  INTRINSIC min
1273  INTEGER :: max1
1274  INTEGER :: max2
1275  INTEGER :: max3
1276  INTEGER :: max4
1277  INTEGER :: max5
1278  INTEGER :: max6
1279  INTEGER :: max7
1280  INTEGER :: max8
1281  INTEGER :: max9
1282  INTEGER :: max10
1283  INTEGER :: max11
1284  INTEGER :: max12
1285  INTEGER :: max13
1286  INTEGER :: max14
1287  INTEGER :: max15
1288  INTEGER :: min1
1289  INTEGER :: min2
1290  INTEGER :: min3
1291  INTEGER :: min4
1292  INTEGER :: min5
1293  INTEGER :: min6
1294  INTEGER :: min7
1295  INTEGER :: min8
1296  INTEGER :: min9
1297  INTEGER :: min10
1298  INTEGER :: min11
1299  INTEGER :: min12
1300  INTEGER :: min13
1301  INTEGER :: min14
1302  INTEGER :: min15
1303  REAL(kind=kind_real) :: result1
1304  REAL(kind=kind_real) :: result1_ad
1305  REAL(kind=kind_real) :: result2
1306  REAL(kind=kind_real) :: result2_ad
1307  REAL(kind=kind_real) :: result3
1308  REAL(kind=kind_real) :: result3_ad
1309  REAL(kind=kind_real) :: temp_ad
1310  REAL(kind=kind_real) :: temp_ad0
1311  REAL(kind=kind_real) :: temp_ad1
1312  REAL(kind=kind_real) :: temp_ad2
1313  REAL(kind=kind_real) :: temp_ad3
1314  real(kind=kind_real) :: temp_ad4
1315  real(kind=kind_real) :: temp_ad5
1316  real(kind=kind_real) :: temp_ad6
1317  real(kind=kind_real) :: temp_ad7
1318  REAL(kind=kind_real) :: temp_ad8
1319  real(kind=kind_real) :: temp_ad9
1320  real(kind=kind_real) :: temp_ad10
1321  real(kind=kind_real) :: temp_ad11
1322  real(kind=kind_real) :: temp_ad12
1323  REAL(kind=kind_real) :: temp_ad13
1324  real(kind=kind_real) :: temp_ad14
1325  real(kind=kind_real) :: temp_ad15
1326  real(kind=kind_real) :: temp_ad16
1327  real(kind=kind_real) :: temp_ad17
1328  REAL(kind=kind_real) :: temp_ad18
1329  real(kind=kind_real) :: temp_ad19
1330  real(kind=kind_real) :: temp_ad20
1331  real(kind=kind_real) :: temp_ad21
1332  real(kind=kind_real) :: temp_ad22
1333  INTEGER :: ad_from
1334  INTEGER :: ad_to
1335  INTEGER :: ad_from0
1336  INTEGER :: ad_to0
1337  INTEGER :: ad_from1
1338  INTEGER :: ad_to1
1339  INTEGER :: ad_from2
1340  INTEGER :: ad_to2
1341  INTEGER :: ad_from3
1342  INTEGER :: ad_to3
1343  INTEGER :: branch
1344  IF (1 .LT. is - 1) THEN
1345  is1 = is - 1
1346  ELSE
1347  is1 = 1
1348  END IF
1349  IF (1 .LT. js - 1) THEN
1350  js1 = js - 1
1351  ELSE
1352  js1 = 1
1353  END IF
1354  IF (2 .LT. is) THEN
1355  is2 = is
1356  ELSE
1357  is2 = 2
1358  END IF
1359  IF (2 .LT. js) THEN
1360  js2 = js
1361  ELSE
1362  js2 = 2
1363  END IF
1364  IF (npx - 1 .GT. ie + 1) THEN
1365  ie1 = ie + 1
1366  ELSE
1367  ie1 = npx - 1
1368  END IF
1369  IF (npy - 1 .GT. je + 1) THEN
1370  je1 = je + 1
1371  ELSE
1372  je1 = npy - 1
1373  END IF
1374 ! Corners:
1375 ! 3-way extrapolation
1376  IF (geom%sw_corner) THEN
1377  CALL pushcontrol1b(0)
1378  ELSE
1379  CALL pushcontrol1b(1)
1380  END IF
1381  IF (geom%se_corner) THEN
1382  CALL pushcontrol1b(0)
1383  ELSE
1384  CALL pushcontrol1b(1)
1385  END IF
1386  IF (geom%ne_corner) THEN
1387  CALL pushcontrol1b(0)
1388  ELSE
1389  CALL pushcontrol1b(1)
1390  END IF
1391  IF (geom%nw_corner) THEN
1392  p0(1:2) = geom%grid(1, npy, 1:2)
1393  CALL pushcontrol1b(1)
1394  ELSE
1395  CALL pushcontrol1b(0)
1396  END IF
1397  IF (1 .LT. js - 2) THEN
1398  max1 = js - 2
1399  ELSE
1400  max1 = 1
1401  END IF
1402  IF (npy - 1 .GT. je + 2) THEN
1403  min1 = je + 2
1404  ELSE
1405  min1 = npy - 1
1406  END IF
1407 !------------
1408 ! X-Interior:
1409 !------------
1410  DO j=max1,min1
1411  IF (3 .LT. is) THEN
1412  max2 = is
1413  ELSE
1414  max2 = 3
1415  END IF
1416  IF (npx - 2 .GT. ie + 1) THEN
1417  min2 = ie + 1
1418  ELSE
1419  min2 = npx - 2
1420  END IF
1421  ad_from = max2
1422  i = min2 + 1
1423  CALL pushinteger4(i - 1)
1424  CALL pushinteger4(ad_from)
1425  END DO
1426 ! *** West Edges:
1427  IF (is .EQ. 1) THEN
1428  IF (1 .LT. js - 2) THEN
1429  max3 = js - 2
1430  ELSE
1431  max3 = 1
1432  END IF
1433  IF (npy - 1 .GT. je + 2) THEN
1434  min3 = je + 2
1435  ELSE
1436  min3 = npy - 1
1437  END IF
1438  CALL pushcontrol1b(0)
1439  ELSE
1440  CALL pushcontrol1b(1)
1441  END IF
1442 ! East Edges:
1443  IF (ie + 1 .EQ. npx) THEN
1444  IF (1 .LT. js - 2) THEN
1445  max4 = js - 2
1446  ELSE
1447  max4 = 1
1448  END IF
1449  IF (npy - 1 .GT. je + 2) THEN
1450  min4 = je + 2
1451  ELSE
1452  min4 = npy - 1
1453  END IF
1454  CALL pushcontrol1b(1)
1455  ELSE
1456  CALL pushcontrol1b(0)
1457  END IF
1458  IF (3 .LT. js) THEN
1459  max5 = js
1460  ELSE
1461  max5 = 3
1462  END IF
1463  IF (npy - 2 .GT. je + 1) THEN
1464  min5 = je + 1
1465  ELSE
1466  min5 = npy - 2
1467  END IF
1468 !------------
1469 ! Y-Interior:
1470 !------------
1471  DO j=max5,min5
1472  IF (1 .LT. is - 2) THEN
1473  max6 = is - 2
1474  ELSE
1475  max6 = 1
1476  END IF
1477  IF (npx - 1 .GT. ie + 2) THEN
1478  min6 = ie + 2
1479  ELSE
1480  min6 = npx - 1
1481  END IF
1482  ad_from0 = max6
1483  i = min6 + 1
1484  CALL pushinteger4(i - 1)
1485  CALL pushinteger4(ad_from0)
1486  END DO
1487 ! South Edges:
1488  IF (js .EQ. 1) THEN
1489  IF (1 .LT. is - 2) THEN
1490  max7 = is - 2
1491  ELSE
1492  max7 = 1
1493  END IF
1494  IF (npx - 1 .GT. ie + 2) THEN
1495  min7 = ie + 2
1496  ELSE
1497  min7 = npx - 1
1498  END IF
1499  CALL pushcontrol1b(0)
1500  ELSE
1501  CALL pushcontrol1b(1)
1502  END IF
1503 ! North Edges:
1504  IF (je + 1 .EQ. npy) THEN
1505  IF (1 .LT. is - 2) THEN
1506  max8 = is - 2
1507  ELSE
1508  max8 = 1
1509  END IF
1510  IF (npx - 1 .GT. ie + 2) THEN
1511  min8 = ie + 2
1512  ELSE
1513  min8 = npx - 1
1514  END IF
1515  CALL pushcontrol1b(1)
1516  ELSE
1517  CALL pushcontrol1b(0)
1518  END IF
1519  IF (3 .LT. js) THEN
1520  max9 = js
1521  ELSE
1522  max9 = 3
1523  END IF
1524  IF (npy - 2 .GT. je + 1) THEN
1525  min9 = je + 1
1526  ELSE
1527  min9 = npy - 2
1528  END IF
1529 !--------------------------------------
1530  DO j=max9,min9
1531  IF (2 .LT. is) THEN
1532  max10 = is
1533  ELSE
1534  max10 = 2
1535  END IF
1536  IF (npx - 1 .GT. ie + 1) THEN
1537  min10 = ie + 1
1538  ELSE
1539  min10 = npx - 1
1540  END IF
1541  ad_from1 = max10
1542  i = min10 + 1
1543  CALL pushinteger4(i - 1)
1544  CALL pushinteger4(ad_from1)
1545  END DO
1546  IF (js .EQ. 1) THEN
1547  IF (2 .LT. is) THEN
1548  max11 = is
1549  ELSE
1550  max11 = 2
1551  END IF
1552  IF (npx - 1 .GT. ie + 1) THEN
1553  min11 = ie + 1
1554  ELSE
1555  min11 = npx - 1
1556  END IF
1557  CALL pushcontrol1b(0)
1558  ELSE
1559  CALL pushcontrol1b(1)
1560  END IF
1561  IF (je + 1 .EQ. npy) THEN
1562  IF (2 .LT. is) THEN
1563  max12 = is
1564  ELSE
1565  max12 = 2
1566  END IF
1567  IF (npx - 1 .GT. ie + 1) THEN
1568  min12 = ie + 1
1569  ELSE
1570  min12 = npx - 1
1571  END IF
1572  CALL pushcontrol1b(1)
1573  ELSE
1574  CALL pushcontrol1b(0)
1575  END IF
1576  IF (2 .LT. js) THEN
1577  max13 = js
1578  ELSE
1579  max13 = 2
1580  END IF
1581  IF (npy - 1 .GT. je + 1) THEN
1582  min13 = je + 1
1583  ELSE
1584  min13 = npy - 1
1585  END IF
1586  DO j=max13,min13
1587  IF (3 .LT. is) THEN
1588  max14 = is
1589  ELSE
1590  max14 = 3
1591  END IF
1592  IF (npx - 2 .GT. ie + 1) THEN
1593  min14 = ie + 1
1594  ELSE
1595  min14 = npx - 2
1596  END IF
1597  ad_from2 = max14
1598  i = min14 + 1
1599  CALL pushinteger4(i - 1)
1600  CALL pushinteger4(ad_from2)
1601  IF (is .EQ. 1) THEN
1602  CALL pushcontrol1b(0)
1603  ELSE
1604  CALL pushcontrol1b(1)
1605  END IF
1606  IF (ie + 1 .EQ. npx) THEN
1607  CALL pushcontrol1b(1)
1608  ELSE
1609  CALL pushcontrol1b(0)
1610  END IF
1611  IF (2 .LT. is) THEN
1612  max15 = is
1613  ELSE
1614  max15 = 2
1615  END IF
1616  IF (npx - 1 .GT. ie + 1) THEN
1617  min15 = ie + 1
1618  ELSE
1619  min15 = npx - 1
1620  END IF
1621  ad_from3 = max15
1622  i = min15 + 1
1623  CALL pushinteger4(i - 1)
1624  CALL pushinteger4(ad_from3)
1625  END DO
1626  qy_ad = 0.0_8
1627  qxx_ad = 0.0_8
1628  qyy_ad = 0.0_8
1629  DO j=min13,max13,-1
1630  CALL popinteger4(ad_from3)
1631  CALL popinteger4(ad_to3)
1632  DO i=ad_to3,ad_from3,-1
1633  qxx_ad(i, j) = qxx_ad(i, j) + 0.5*qout_ad(i, j)
1634  qyy_ad(i, j) = qyy_ad(i, j) + 0.5*qout_ad(i, j)
1635  qout_ad(i, j) = 0.0_8
1636  END DO
1637  CALL popcontrol1b(branch)
1638  IF (branch .NE. 0) THEN
1639  qy_ad(npx-2, j) = qy_ad(npx-2, j) + c1*qyy_ad(npx-1, j)
1640  qy_ad(npx-1, j) = qy_ad(npx-1, j) + c1*qyy_ad(npx-1, j)
1641  qout_ad(npx, j) = qout_ad(npx, j) + c2*qyy_ad(npx-1, j)
1642  qyy_ad(npx-2, j) = qyy_ad(npx-2, j) + c2*qyy_ad(npx-1, j)
1643  qyy_ad(npx-1, j) = 0.0_8
1644  END IF
1645  CALL popcontrol1b(branch)
1646  IF (branch .EQ. 0) THEN
1647  qy_ad(1, j) = qy_ad(1, j) + c1*qyy_ad(2, j)
1648  qy_ad(2, j) = qy_ad(2, j) + c1*qyy_ad(2, j)
1649  qout_ad(1, j) = qout_ad(1, j) + c2*qyy_ad(2, j)
1650  qyy_ad(3, j) = qyy_ad(3, j) + c2*qyy_ad(2, j)
1651  qyy_ad(2, j) = 0.0_8
1652  END IF
1653  CALL popinteger4(ad_from2)
1654  CALL popinteger4(ad_to2)
1655  DO i=ad_to2,ad_from2,-1
1656  qy_ad(i-2, j) = qy_ad(i-2, j) + a2*qyy_ad(i, j)
1657  qy_ad(i+1, j) = qy_ad(i+1, j) + a2*qyy_ad(i, j)
1658  qy_ad(i-1, j) = qy_ad(i-1, j) + a1*qyy_ad(i, j)
1659  qy_ad(i, j) = qy_ad(i, j) + a1*qyy_ad(i, j)
1660  qyy_ad(i, j) = 0.0_8
1661  END DO
1662  END DO
1663  CALL popcontrol1b(branch)
1664  IF (branch .EQ. 0) THEN
1665  qx_ad = 0.0_8
1666  ELSE
1667  qx_ad = 0.0_8
1668  DO i=min12,max12,-1
1669  qx_ad(i, npy-2) = qx_ad(i, npy-2) + c1*qxx_ad(i, npy-1)
1670  qx_ad(i, npy-1) = qx_ad(i, npy-1) + c1*qxx_ad(i, npy-1)
1671  qout_ad(i, npy) = qout_ad(i, npy) + c2*qxx_ad(i, npy-1)
1672  qxx_ad(i, npy-2) = qxx_ad(i, npy-2) + c2*qxx_ad(i, npy-1)
1673  qxx_ad(i, npy-1) = 0.0_8
1674  END DO
1675  END IF
1676  CALL popcontrol1b(branch)
1677  IF (branch .EQ. 0) THEN
1678  DO i=min11,max11,-1
1679  qx_ad(i, 1) = qx_ad(i, 1) + c1*qxx_ad(i, 2)
1680  qx_ad(i, 2) = qx_ad(i, 2) + c1*qxx_ad(i, 2)
1681  qout_ad(i, 1) = qout_ad(i, 1) + c2*qxx_ad(i, 2)
1682  qxx_ad(i, 3) = qxx_ad(i, 3) + c2*qxx_ad(i, 2)
1683  qxx_ad(i, 2) = 0.0_8
1684  END DO
1685  END IF
1686  DO j=min9,max9,-1
1687  CALL popinteger4(ad_from1)
1688  CALL popinteger4(ad_to1)
1689  DO i=ad_to1,ad_from1,-1
1690  qx_ad(i, j-2) = qx_ad(i, j-2) + a2*qxx_ad(i, j)
1691  qx_ad(i, j+1) = qx_ad(i, j+1) + a2*qxx_ad(i, j)
1692  qx_ad(i, j-1) = qx_ad(i, j-1) + a1*qxx_ad(i, j)
1693  qx_ad(i, j) = qx_ad(i, j) + a1*qxx_ad(i, j)
1694  qxx_ad(i, j) = 0.0_8
1695  END DO
1696  END DO
1697  CALL popcontrol1b(branch)
1698  IF (branch .EQ. 0) THEN
1699  q1_ad = 0.0_8
1700  ELSE
1701  DO i=min8,max8,-1
1702  g_in = geom%dya(i, npy-2)/geom%dya(i, npy-1)
1703  temp_ad19 = qy_ad(i, npy-1)/(g_in*2.+2.)
1704  qin_ad(i, npy-2) = qin_ad(i, npy-2) + 3.*temp_ad19
1705  qy_ad(i, npy) = qy_ad(i, npy) - g_in*temp_ad19
1706  qy_ad(i, npy-2) = qy_ad(i, npy-2) - temp_ad19
1707  qy_ad(i, npy-1) = 0.0_8
1708  g_ou = geom%dya(i, npy+1)/geom%dya(i, npy)
1709  temp_ad21 = 0.5*qy_ad(i, npy)
1710  temp_ad20 = temp_ad21/(g_in+1.)
1711  qin_ad(i, npy-1) = qin_ad(i, npy-1) + (g_in+2.)*temp_ad20 + 3.*&
1712 & g_in*temp_ad19
1713  temp_ad22 = temp_ad21/(g_ou+1.)
1714  qin_ad(i, npy-2) = qin_ad(i, npy-2) - temp_ad20
1715  qin_ad(i, npy) = qin_ad(i, npy) + (g_ou+2.)*temp_ad22
1716  qin_ad(i, npy+1) = qin_ad(i, npy+1) - temp_ad22
1717  qy_ad(i, npy) = 0.0_8
1718  END DO
1719  q1_ad = 0.0_8
1720  DO i=ie1,is2,-1
1721  q1_ad(i-1) = q1_ad(i-1) + geom%edge_n(i)*qout_ad(i, npy)
1722  q1_ad(i) = q1_ad(i) + (1.-geom%edge_n(i))*qout_ad(i, npy)
1723  qout_ad(i, npy) = 0.0_8
1724  END DO
1725  DO i=ie1,is1,-1
1726  temp_ad18 = q1_ad(i)/(geom%dya(i, npy-1)+geom%dya(i, npy))
1727  qin_ad(i, npy-1) = qin_ad(i, npy-1) + geom%dya(i, npy)*temp_ad18
1728  qin_ad(i, npy) = qin_ad(i, npy) + geom%dya(i, npy-1)*temp_ad18
1729  q1_ad(i) = 0.0_8
1730  END DO
1731  END IF
1732  CALL popcontrol1b(branch)
1733  IF (branch .EQ. 0) THEN
1734  DO i=min7,max7,-1
1735  g_in = geom%dya(i, 2)/geom%dya(i, 1)
1736  temp_ad14 = qy_ad(i, 2)/(g_in*2.+2.)
1737  qin_ad(i, 1) = qin_ad(i, 1) + 3.*g_in*temp_ad14
1738  qin_ad(i, 2) = qin_ad(i, 2) + 3.*temp_ad14
1739  qy_ad(i, 1) = qy_ad(i, 1) - g_in*temp_ad14
1740  qy_ad(i, 3) = qy_ad(i, 3) - temp_ad14
1741  qy_ad(i, 2) = 0.0_8
1742  g_ou = geom%dya(i, -1)/geom%dya(i, 0)
1743  temp_ad15 = 0.5*qy_ad(i, 1)
1744  temp_ad16 = temp_ad15/(g_in+1.)
1745  temp_ad17 = temp_ad15/(g_ou+1.)
1746  qin_ad(i, 1) = qin_ad(i, 1) + (g_in+2.)*temp_ad16
1747  qin_ad(i, 2) = qin_ad(i, 2) - temp_ad16
1748  qin_ad(i, 0) = qin_ad(i, 0) + (g_ou+2.)*temp_ad17
1749  qin_ad(i, -1) = qin_ad(i, -1) - temp_ad17
1750  qy_ad(i, 1) = 0.0_8
1751  END DO
1752  DO i=ie1,is2,-1
1753  q1_ad(i-1) = q1_ad(i-1) + geom%edge_s(i)*qout_ad(i, 1)
1754  q1_ad(i) = q1_ad(i) + (1.-geom%edge_s(i))*qout_ad(i, 1)
1755  qout_ad(i, 1) = 0.0_8
1756  END DO
1757  DO i=ie1,is1,-1
1758  temp_ad13 = q1_ad(i)/(geom%dya(i, 0)+geom%dya(i, 1))
1759  qin_ad(i, 0) = qin_ad(i, 0) + geom%dya(i, 1)*temp_ad13
1760  qin_ad(i, 1) = qin_ad(i, 1) + geom%dya(i, 0)*temp_ad13
1761  q1_ad(i) = 0.0_8
1762  END DO
1763  END IF
1764  DO j=min5,max5,-1
1765  CALL popinteger4(ad_from0)
1766  CALL popinteger4(ad_to0)
1767  DO i=ad_to0,ad_from0,-1
1768  qin_ad(i, j-2) = qin_ad(i, j-2) + b2*qy_ad(i, j)
1769  qin_ad(i, j+1) = qin_ad(i, j+1) + b2*qy_ad(i, j)
1770  qin_ad(i, j-1) = qin_ad(i, j-1) + b1*qy_ad(i, j)
1771  qin_ad(i, j) = qin_ad(i, j) + b1*qy_ad(i, j)
1772  qy_ad(i, j) = 0.0_8
1773  END DO
1774  END DO
1775  CALL popcontrol1b(branch)
1776  IF (branch .EQ. 0) THEN
1777  q2_ad = 0.0_8
1778  ELSE
1779  DO j=min4,max4,-1
1780  g_in = geom%dxa(npx-2, j)/geom%dxa(npx-1, j)
1781  temp_ad9 = qx_ad(npx-1, j)/(g_in*2.+2.)
1782  qin_ad(npx-2, j) = qin_ad(npx-2, j) + 3.*temp_ad9
1783  qx_ad(npx, j) = qx_ad(npx, j) - g_in*temp_ad9
1784  qx_ad(npx-2, j) = qx_ad(npx-2, j) - temp_ad9
1785  qx_ad(npx-1, j) = 0.0_8
1786  g_ou = geom%dxa(npx+1, j)/geom%dxa(npx, j)
1787  temp_ad11 = 0.5*qx_ad(npx, j)
1788  temp_ad10 = temp_ad11/(g_in+1.)
1789  qin_ad(npx-1, j) = qin_ad(npx-1, j) + (g_in+2.)*temp_ad10 + 3.*&
1790 & g_in*temp_ad9
1791  temp_ad12 = temp_ad11/(g_ou+1.)
1792  qin_ad(npx-2, j) = qin_ad(npx-2, j) - temp_ad10
1793  qin_ad(npx, j) = qin_ad(npx, j) + (g_ou+2.)*temp_ad12
1794  qin_ad(npx+1, j) = qin_ad(npx+1, j) - temp_ad12
1795  qx_ad(npx, j) = 0.0_8
1796  END DO
1797  q2_ad = 0.0_8
1798  DO j=je1,js2,-1
1799  q2_ad(j-1) = q2_ad(j-1) + geom%edge_e(j)*qout_ad(npx, j)
1800  q2_ad(j) = q2_ad(j) + (1.-geom%edge_e(j))*qout_ad(npx, j)
1801  qout_ad(npx, j) = 0.0_8
1802  END DO
1803  DO j=je1,js1,-1
1804  temp_ad8 = q2_ad(j)/(geom%dxa(npx-1, j)+geom%dxa(npx, j))
1805  qin_ad(npx-1, j) = qin_ad(npx-1, j) + geom%dxa(npx, j)*temp_ad8
1806  qin_ad(npx, j) = qin_ad(npx, j) + geom%dxa(npx-1, j)*temp_ad8
1807  q2_ad(j) = 0.0_8
1808  END DO
1809  END IF
1810  CALL popcontrol1b(branch)
1811  IF (branch .EQ. 0) THEN
1812  DO j=min3,max3,-1
1813  g_in = geom%dxa(2, j)/geom%dxa(1, j)
1814  temp_ad4 = qx_ad(2, j)/(g_in*2.+2.)
1815  qin_ad(1, j) = qin_ad(1, j) + 3.*g_in*temp_ad4
1816  qin_ad(2, j) = qin_ad(2, j) + 3.*temp_ad4
1817  qx_ad(1, j) = qx_ad(1, j) - g_in*temp_ad4
1818  qx_ad(3, j) = qx_ad(3, j) - temp_ad4
1819  qx_ad(2, j) = 0.0_8
1820  g_ou = geom%dxa(-1, j)/geom%dxa(0, j)
1821  temp_ad5 = 0.5*qx_ad(1, j)
1822  temp_ad6 = temp_ad5/(g_in+1.)
1823  temp_ad7 = temp_ad5/(g_ou+1.)
1824  qin_ad(1, j) = qin_ad(1, j) + (g_in+2.)*temp_ad6
1825  qin_ad(2, j) = qin_ad(2, j) - temp_ad6
1826  qin_ad(0, j) = qin_ad(0, j) + (g_ou+2.)*temp_ad7
1827  qin_ad(-1, j) = qin_ad(-1, j) - temp_ad7
1828  qx_ad(1, j) = 0.0_8
1829  END DO
1830  DO j=je1,js2,-1
1831  q2_ad(j-1) = q2_ad(j-1) + geom%edge_w(j)*qout_ad(1, j)
1832  q2_ad(j) = q2_ad(j) + (1.-geom%edge_w(j))*qout_ad(1, j)
1833  qout_ad(1, j) = 0.0_8
1834  END DO
1835  DO j=je1,js1,-1
1836  temp_ad3 = q2_ad(j)/(geom%dxa(0, j)+geom%dxa(1, j))
1837  qin_ad(0, j) = qin_ad(0, j) + geom%dxa(1, j)*temp_ad3
1838  qin_ad(1, j) = qin_ad(1, j) + geom%dxa(0, j)*temp_ad3
1839  q2_ad(j) = 0.0_8
1840  END DO
1841  END IF
1842  DO j=min1,max1,-1
1843  CALL popinteger4(ad_from)
1844  CALL popinteger4(ad_to)
1845  DO i=ad_to,ad_from,-1
1846  qin_ad(i-2, j) = qin_ad(i-2, j) + b2*qx_ad(i, j)
1847  qin_ad(i+1, j) = qin_ad(i+1, j) + b2*qx_ad(i, j)
1848  qin_ad(i-1, j) = qin_ad(i-1, j) + b1*qx_ad(i, j)
1849  qin_ad(i, j) = qin_ad(i, j) + b1*qx_ad(i, j)
1850  qx_ad(i, j) = 0.0_8
1851  END DO
1852  END DO
1853  CALL popcontrol1b(branch)
1854  IF (branch .NE. 0) THEN
1855  temp_ad2 = r3*qout_ad(1, npy)
1856  result1_ad = temp_ad2
1857  result2_ad = temp_ad2
1858  result3_ad = temp_ad2
1859  qout_ad(1, npy) = 0.0_8
1860  CALL extrap_corner_adm(p0, geom%agrid(1, npy, 1:2), geom%agrid(2, &
1861 & npy+1, 1:2), qin(1, npy), qin_ad(1, npy), qin(2, &
1862 & npy+1), qin_ad(2, npy+1), result3_ad)
1863  CALL extrap_corner_adm(p0, geom%agrid(0, npy-1, 1:2), geom%agrid(-&
1864 & 1, npy-2, 1:2), qin(0, npy-1), qin_ad(0, npy-1), &
1865 & qin(-1, npy-2), qin_ad(-1, npy-2), result2_ad)
1866  CALL extrap_corner_adm(p0, geom%agrid(1, npy-1, 1:2), geom%agrid(2&
1867 & , npy-2, 1:2), qin(1, npy-1), qin_ad(1, npy-1), &
1868 & qin(2, npy-2), qin_ad(2, npy-2), result1_ad)
1869  END IF
1870  CALL popcontrol1b(branch)
1871  IF (branch .EQ. 0) THEN
1872  temp_ad1 = r3*qout_ad(npx, npy)
1873  result1_ad = temp_ad1
1874  result2_ad = temp_ad1
1875  result3_ad = temp_ad1
1876  qout_ad(npx, npy) = 0.0_8
1877  p0(1:2) = geom%grid(npx, npy, 1:2)
1878  CALL extrap_corner_adm(p0, geom%agrid(npx-1, npy, 1:2), geom%agrid&
1879 & (npx-2, npy+1, 1:2), qin(npx-1, npy), qin_ad(npx-&
1880 & 1, npy), qin(npx-2, npy+1), qin_ad(npx-2, npy+1)&
1881 & , result3_ad)
1882  CALL extrap_corner_adm(p0, geom%agrid(npx, npy-1, 1:2), geom%agrid&
1883 & (npx+1, npy-2, 1:2), qin(npx, npy-1), qin_ad(npx&
1884 & , npy-1), qin(npx+1, npy-2), qin_ad(npx+1, npy-2)&
1885 & , result2_ad)
1886  CALL extrap_corner_adm(p0, geom%agrid(npx-1, npy-1, 1:2), geom%&
1887 & agrid(npx-2, npy-2, 1:2), qin(npx-1, npy-1), &
1888 & qin_ad(npx-1, npy-1), qin(npx-2, npy-2), qin_ad(&
1889 & npx-2, npy-2), result1_ad)
1890  END IF
1891  CALL popcontrol1b(branch)
1892  IF (branch .EQ. 0) THEN
1893  temp_ad0 = r3*qout_ad(npx, 1)
1894  result1_ad = temp_ad0
1895  result2_ad = temp_ad0
1896  result3_ad = temp_ad0
1897  qout_ad(npx, 1) = 0.0_8
1898  p0(1:2) = geom%grid(npx, 1, 1:2)
1899  CALL extrap_corner_adm(p0, geom%agrid(npx, 1, 1:2), geom%agrid(npx&
1900 & +1, 2, 1:2), qin(npx, 1), qin_ad(npx, 1), qin(npx&
1901 & +1, 2), qin_ad(npx+1, 2), result3_ad)
1902  CALL extrap_corner_adm(p0, geom%agrid(npx-1, 0, 1:2), geom%agrid(&
1903 & npx-2, -1, 1:2), qin(npx-1, 0), qin_ad(npx-1, 0)&
1904 & , qin(npx-2, -1), qin_ad(npx-2, -1), result2_ad)
1905  CALL extrap_corner_adm(p0, geom%agrid(npx-1, 1, 1:2), geom%agrid(&
1906 & npx-2, 2, 1:2), qin(npx-1, 1), qin_ad(npx-1, 1), &
1907 & qin(npx-2, 2), qin_ad(npx-2, 2), result1_ad)
1908  END IF
1909  CALL popcontrol1b(branch)
1910  IF (branch .EQ. 0) THEN
1911  temp_ad = r3*qout_ad(1, 1)
1912  result1_ad = temp_ad
1913  result2_ad = temp_ad
1914  result3_ad = temp_ad
1915  p0(1:2) = geom%grid(1, 1, 1:2)
1916  CALL extrap_corner_adm(p0, geom%agrid(1, 0, 1:2), geom%agrid(2, -1&
1917 & , 1:2), qin(1, 0), qin_ad(1, 0), qin(2, -1), &
1918 & qin_ad(2, -1), result3_ad)
1919  CALL extrap_corner_adm(p0, geom%agrid(0, 1, 1:2), geom%agrid(-1, 2&
1920 & , 1:2), qin(0, 1), qin_ad(0, 1), qin(-1, 2), &
1921 & qin_ad(-1, 2), result2_ad)
1922  CALL extrap_corner_adm(p0, geom%agrid(1, 1, 1:2), geom%agrid(2, 2&
1923 & , 1:2), qin(1, 1), qin_ad(1, 1), qin(2, 2), &
1924 & qin_ad(2, 2), result1_ad)
1925  END IF
1926  END SUBROUTINE a2b_ord4_adm
1927 
1928 !----------------------------------------------------------------------------
1929 
1930  SUBROUTINE extrap_corner_adm(p0, p1, p2, q1, q1_ad, q2, q2_ad, &
1931 & extrap_corner_ad)
1932  IMPLICIT NONE
1933  REAL(kind=kind_real), DIMENSION(2), INTENT(IN) :: p0, p1, p2
1934  REAL(kind=kind_real), INTENT(IN) :: q1, q2
1935  REAL(kind=kind_real) :: q1_ad, q2_ad
1936  REAL(kind=kind_real) :: x1, x2
1937  REAL(kind=kind_real) :: temp_ad
1938  REAL(kind=kind_real) :: extrap_corner_ad
1939  REAL(kind=kind_real) :: extrap_corner
1940  x1 = great_circle_dist(p1, p0)
1941  x2 = great_circle_dist(p2, p0)
1942  temp_ad = x1*extrap_corner_ad/(x2-x1)
1943  q1_ad = q1_ad + temp_ad + extrap_corner_ad
1944  q2_ad = q2_ad - temp_ad
1945  END SUBROUTINE extrap_corner_adm
1946 
1947 ! ------------------------------------------------------------------------------
1948 
1949 subroutine a2d(geom, ua, va, ud, vd)
1950 
1951 use mpp_domains_mod, only: mpp_update_domains
1952 
1953  implicit none
1954 
1955  type(fv3jedi_geom), intent(inout) :: geom
1956  real(kind=kind_real), intent(in) :: ua(geom%isc:geom%iec,geom%jsc:geom%jec,geom%npz)
1957  real(kind=kind_real), intent(in) :: va(geom%isc:geom%iec,geom%jsc:geom%jec,geom%npz)
1958  real(kind=kind_real), intent(inout) :: ud(geom%isc:geom%iec ,geom%jsc:geom%jec+1,geom%npz)
1959  real(kind=kind_real), intent(inout) :: vd(geom%isc:geom%iec+1,geom%jsc:geom%jec ,geom%npz)
1960 
1961  integer :: is ,ie , js ,je
1962  integer :: npx, npy, npz
1963  integer :: i,j,k, im2,jm2
1964 
1965  real(kind=kind_real) :: uatemp(geom%isd:geom%ied,geom%jsd:geom%jed,geom%npz)
1966  real(kind=kind_real) :: vatemp(geom%isd:geom%ied,geom%jsd:geom%jed,geom%npz)
1967 
1968  real(kind=kind_real) :: v3(geom%isc-1:geom%iec+1,geom%jsc-1:geom%jec+1,3)
1969  real(kind=kind_real) :: ue(geom%isc-1:geom%iec+1,geom%jsc :geom%jec+1,3) ! 3D winds at edges
1970  real(kind=kind_real) :: ve(geom%isc :geom%iec+1,geom%jsc-1:geom%jec+1,3) ! 3D winds at edges
1971  real(kind=kind_real), dimension(geom%isc:geom%iec):: ut1, ut2, ut3
1972  real(kind=kind_real), dimension(geom%jsc:geom%jec):: vt1, vt2, vt3
1973 
1974  npx = geom%npx
1975  npy = geom%npy
1976  npz = geom%npz
1977  is = geom%isc
1978  ie = geom%iec
1979  js = geom%jsc
1980  je = geom%jec
1981 
1982  im2 = (npx-1)/2
1983  jm2 = (npy-1)/2
1984 
1985  uatemp(:,:,:) = 0.0
1986  vatemp(:,:,:) = 0.0
1987 
1988  uatemp(is:ie,js:je,:) = ua
1989  vatemp(is:ie,js:je,:) = va
1990 
1991  call mpp_update_domains(uatemp, geom%domain, complete=.true.)
1992  call mpp_update_domains(vatemp, geom%domain, complete=.true.)
1993 
1994  do k=1, npz
1995 
1996  do j=js-1,je+1
1997  do i=is-1,ie+1
1998  v3(i,j,1) = uatemp(i,j,k)*geom%vlon(i,j,1) + vatemp(i,j,k)*geom%vlat(i,j,1)
1999  v3(i,j,2) = uatemp(i,j,k)*geom%vlon(i,j,2) + vatemp(i,j,k)*geom%vlat(i,j,2)
2000  v3(i,j,3) = uatemp(i,j,k)*geom%vlon(i,j,3) + vatemp(i,j,k)*geom%vlat(i,j,3)
2001  enddo
2002  enddo
2003 
2004  do j=js,je+1
2005  do i=is-1,ie+1
2006  ue(i,j,1) = v3(i,j-1,1) + v3(i,j,1)
2007  ue(i,j,2) = v3(i,j-1,2) + v3(i,j,2)
2008  ue(i,j,3) = v3(i,j-1,3) + v3(i,j,3)
2009  enddo
2010  enddo
2011 
2012  do j=js-1,je+1
2013  do i=is,ie+1
2014  ve(i,j,1) = v3(i-1,j,1) + v3(i,j,1)
2015  ve(i,j,2) = v3(i-1,j,2) + v3(i,j,2)
2016  ve(i,j,3) = v3(i-1,j,3) + v3(i,j,3)
2017  enddo
2018  enddo
2019 
2020  if ( is==1 ) then
2021  i = 1
2022  do j=js,je
2023  if ( j>jm2 ) then
2024  vt1(j) = geom%edge_vect_w(j)*ve(i,j-1,1)+(1.-geom%edge_vect_w(j))*ve(i,j,1)
2025  vt2(j) = geom%edge_vect_w(j)*ve(i,j-1,2)+(1.-geom%edge_vect_w(j))*ve(i,j,2)
2026  vt3(j) = geom%edge_vect_w(j)*ve(i,j-1,3)+(1.-geom%edge_vect_w(j))*ve(i,j,3)
2027  else
2028  vt1(j) = geom%edge_vect_w(j)*ve(i,j+1,1)+(1.-geom%edge_vect_w(j))*ve(i,j,1)
2029  vt2(j) = geom%edge_vect_w(j)*ve(i,j+1,2)+(1.-geom%edge_vect_w(j))*ve(i,j,2)
2030  vt3(j) = geom%edge_vect_w(j)*ve(i,j+1,3)+(1.-geom%edge_vect_w(j))*ve(i,j,3)
2031  endif
2032  enddo
2033  do j=js,je
2034  ve(i,j,1) = vt1(j)
2035  ve(i,j,2) = vt2(j)
2036  ve(i,j,3) = vt3(j)
2037  enddo
2038  endif
2039 
2040  if ( (ie+1)==npx ) then
2041  i = npx
2042  do j=js,je
2043  if ( j>jm2 ) then
2044  vt1(j) = geom%edge_vect_e(j)*ve(i,j-1,1)+(1.-geom%edge_vect_e(j))*ve(i,j,1)
2045  vt2(j) = geom%edge_vect_e(j)*ve(i,j-1,2)+(1.-geom%edge_vect_e(j))*ve(i,j,2)
2046  vt3(j) = geom%edge_vect_e(j)*ve(i,j-1,3)+(1.-geom%edge_vect_e(j))*ve(i,j,3)
2047  else
2048  vt1(j) = geom%edge_vect_e(j)*ve(i,j+1,1)+(1.-geom%edge_vect_e(j))*ve(i,j,1)
2049  vt2(j) = geom%edge_vect_e(j)*ve(i,j+1,2)+(1.-geom%edge_vect_e(j))*ve(i,j,2)
2050  vt3(j) = geom%edge_vect_e(j)*ve(i,j+1,3)+(1.-geom%edge_vect_e(j))*ve(i,j,3)
2051  endif
2052  enddo
2053  do j=js,je
2054  ve(i,j,1) = vt1(j)
2055  ve(i,j,2) = vt2(j)
2056  ve(i,j,3) = vt3(j)
2057  enddo
2058  endif
2059 
2060  if ( js==1 ) then
2061  j = 1
2062  do i=is,ie
2063  if ( i>im2 ) then
2064  ut1(i) = geom%edge_vect_s(i)*ue(i-1,j,1)+(1.-geom%edge_vect_s(i))*ue(i,j,1)
2065  ut2(i) = geom%edge_vect_s(i)*ue(i-1,j,2)+(1.-geom%edge_vect_s(i))*ue(i,j,2)
2066  ut3(i) = geom%edge_vect_s(i)*ue(i-1,j,3)+(1.-geom%edge_vect_s(i))*ue(i,j,3)
2067  else
2068  ut1(i) = geom%edge_vect_s(i)*ue(i+1,j,1)+(1.-geom%edge_vect_s(i))*ue(i,j,1)
2069  ut2(i) = geom%edge_vect_s(i)*ue(i+1,j,2)+(1.-geom%edge_vect_s(i))*ue(i,j,2)
2070  ut3(i) = geom%edge_vect_s(i)*ue(i+1,j,3)+(1.-geom%edge_vect_s(i))*ue(i,j,3)
2071  endif
2072  enddo
2073  do i=is,ie
2074  ue(i,j,1) = ut1(i)
2075  ue(i,j,2) = ut2(i)
2076  ue(i,j,3) = ut3(i)
2077  enddo
2078  endif
2079 
2080  if ( (je+1)==npy ) then
2081  j = npy
2082  do i=is,ie
2083  if ( i>im2 ) then
2084  ut1(i) = geom%edge_vect_n(i)*ue(i-1,j,1)+(1.-geom%edge_vect_n(i))*ue(i,j,1)
2085  ut2(i) = geom%edge_vect_n(i)*ue(i-1,j,2)+(1.-geom%edge_vect_n(i))*ue(i,j,2)
2086  ut3(i) = geom%edge_vect_n(i)*ue(i-1,j,3)+(1.-geom%edge_vect_n(i))*ue(i,j,3)
2087  else
2088  ut1(i) = geom%edge_vect_n(i)*ue(i+1,j,1)+(1.-geom%edge_vect_n(i))*ue(i,j,1)
2089  ut2(i) = geom%edge_vect_n(i)*ue(i+1,j,2)+(1.-geom%edge_vect_n(i))*ue(i,j,2)
2090  ut3(i) = geom%edge_vect_n(i)*ue(i+1,j,3)+(1.-geom%edge_vect_n(i))*ue(i,j,3)
2091  endif
2092  enddo
2093  do i=is,ie
2094  ue(i,j,1) = ut1(i)
2095  ue(i,j,2) = ut2(i)
2096  ue(i,j,3) = ut3(i)
2097  enddo
2098  endif
2099 
2100  do j=js,je+1
2101  do i=is,ie
2102  ud(i,j,k) = 0.5*( ue(i,j,1)*geom%es(1,i,j,1) + &
2103  ue(i,j,2)*geom%es(2,i,j,1) + &
2104  ue(i,j,3)*geom%es(3,i,j,1) )
2105  enddo
2106  enddo
2107 
2108  do j=js,je
2109  do i=is,ie+1
2110  vd(i,j,k) = 0.5*( ve(i,j,1)*geom%ew(1,i,j,2) + &
2111  ve(i,j,2)*geom%ew(2,i,j,2) + &
2112  ve(i,j,3)*geom%ew(3,i,j,2) )
2113  enddo
2114  enddo
2115 
2116  enddo
2117 
2118 end subroutine a2d
2119 
2120 ! ------------------------------------------------------------------------------
2121 
2122 subroutine a2d_ad(geom, ua_ad, va_ad, ud_ad, vd_ad)
2123 
2124 use fv_mp_adm_mod, only: mpp_update_domains_adm
2125 
2126 implicit none
2127 type(fv3jedi_geom), intent(inout) :: geom
2128 real(kind=kind_real), intent(inout) :: ua_ad(geom%isc:geom%iec, geom%jsc:geom%jec, geom%npz)
2129 real(kind=kind_real), intent(inout) :: va_ad(geom%isc:geom%iec, geom%jsc:geom%jec, geom%npz)
2130 real(kind=kind_real), intent(inout) :: ud_ad(geom%isc:geom%iec , geom%jsc:geom%jec+1, geom%npz)
2131 real(kind=kind_real), intent(inout) :: vd_ad(geom%isc:geom%iec+1, geom%jsc:geom%jec , geom%npz)
2132 
2133 integer :: is, ie, js, je
2134 integer :: npx, npy, npz
2135 integer :: i, j, k, im2, jm2
2136 
2137 real(kind=kind_real) :: uatemp(geom%isd:geom%ied, geom%jsd:geom%jed, geom%npz)
2138 real(kind=kind_real) :: vatemp(geom%isd:geom%ied, geom%jsd:geom%jed, geom%npz)
2139 
2140 real(kind=kind_real) :: uatemp_ad(geom%isd:geom%ied, geom%jsd:geom%jed, geom%npz)
2141 real(kind=kind_real) :: vatemp_ad(geom%isd:geom%ied, geom%jsd:geom%jed, geom%npz)
2142 real(kind=kind_real) :: v3_ad(geom%isc-1:geom%iec+1, geom%jsc-1:geom%jec+1, 3)
2143 real(kind=kind_real) :: ue_ad(geom%isc-1:geom%iec+1, geom%jsc:geom%jec+1, 3)
2144 real(kind=kind_real) :: ve_ad(geom%isc:geom%iec+1, geom%jsc-1:geom%jec+1, 3)
2145 real(kind=kind_real), dimension(geom%isc:geom%iec) :: ut1_ad, ut2_ad, ut3_ad
2146 real(kind=kind_real), dimension(geom%jsc:geom%jec) :: vt1_ad, vt2_ad, vt3_ad
2147 real(kind=kind_real) :: temp_ad
2148 real(kind=kind_real) :: temp_ad0
2149 
2150  npx = geom%npx
2151  npy = geom%npy
2152  npz = geom%npz
2153  is = geom%isc
2154  ie = geom%iec
2155  js = geom%jsc
2156  je = geom%jec
2157  im2 = (npx-1)/2
2158  jm2 = (npy-1)/2
2159 
2160  uatemp(:, :, :) = 0.0
2161  vatemp(:, :, :) = 0.0
2162 
2163  v3_ad = 0.0_8
2164  uatemp_ad = 0.0_8
2165  ue_ad = 0.0_8
2166  vt1_ad = 0.0_8
2167  vt2_ad = 0.0_8
2168  vt3_ad = 0.0_8
2169  ve_ad = 0.0_8
2170  vatemp_ad = 0.0_8
2171  ut1_ad = 0.0_8
2172  ut2_ad = 0.0_8
2173  ut3_ad = 0.0_8
2174 
2175  do k=npz,1,-1
2176 
2177  do j=je,js,-1
2178  do i=ie+1,is,-1
2179  temp_ad0 = 0.5*vd_ad(i, j, k)
2180  ve_ad(i, j, 1) = ve_ad(i, j, 1) + geom%ew(1, i, j, 2)*temp_ad0
2181  ve_ad(i, j, 2) = ve_ad(i, j, 2) + geom%ew(2, i, j, 2)*temp_ad0
2182  ve_ad(i, j, 3) = ve_ad(i, j, 3) + geom%ew(3, i, j, 2)*temp_ad0
2183  vd_ad(i, j, k) = 0.0_8
2184  end do
2185  end do
2186 
2187  do j=je+1,js,-1
2188  do i=ie,is,-1
2189  temp_ad = 0.5*ud_ad(i, j, k)
2190  ue_ad(i, j, 1) = ue_ad(i, j, 1) + geom%es(1, i, j, 1)*temp_ad
2191  ue_ad(i, j, 2) = ue_ad(i, j, 2) + geom%es(2, i, j, 1)*temp_ad
2192  ue_ad(i, j, 3) = ue_ad(i, j, 3) + geom%es(3, i, j, 1)*temp_ad
2193  ud_ad(i, j, k) = 0.0_8
2194  end do
2195  end do
2196 
2197  if (je + 1 .eq. npy) then
2198  do i=ie,is,-1
2199  ut3_ad(i) = ut3_ad(i) + ue_ad(i, npy, 3)
2200  ue_ad(i, npy, 3) = 0.0_8
2201  ut2_ad(i) = ut2_ad(i) + ue_ad(i, npy, 2)
2202  ue_ad(i, npy, 2) = 0.0_8
2203  ut1_ad(i) = ut1_ad(i) + ue_ad(i, npy, 1)
2204  ue_ad(i, npy, 1) = 0.0_8
2205  end do
2206  do i=ie,is,-1
2207  if (i .le. im2) then
2208  ue_ad(i+1, npy, 3) = ue_ad(i+1, npy, 3) + geom%edge_vect_n(i)*ut3_ad(i)
2209  ue_ad(i, npy, 3) = ue_ad(i, npy, 3) + (1.-geom%edge_vect_n(i))*ut3_ad(i)
2210  ut3_ad(i) = 0.0_8
2211  ue_ad(i+1, npy, 2) = ue_ad(i+1, npy, 2) + geom%edge_vect_n(i)*ut2_ad(i)
2212  ue_ad(i, npy, 2) = ue_ad(i, npy, 2) + (1.-geom%edge_vect_n(i))*ut2_ad(i)
2213  ut2_ad(i) = 0.0_8
2214  ue_ad(i+1, npy, 1) = ue_ad(i+1, npy, 1) + geom%edge_vect_n(i)*ut1_ad(i)
2215  ue_ad(i, npy, 1) = ue_ad(i, npy, 1) + (1.-geom%edge_vect_n(i))*ut1_ad(i)
2216  ut1_ad(i) = 0.0_8
2217  else
2218  ue_ad(i-1, npy, 3) = ue_ad(i-1, npy, 3) + geom%edge_vect_n(i)*ut3_ad(i)
2219  ue_ad(i, npy, 3) = ue_ad(i, npy, 3) + (1.-geom%edge_vect_n(i))*ut3_ad(i)
2220  ut3_ad(i) = 0.0_8
2221  ue_ad(i-1, npy, 2) = ue_ad(i-1, npy, 2) + geom%edge_vect_n(i)*ut2_ad(i)
2222  ue_ad(i, npy, 2) = ue_ad(i, npy, 2) + (1.-geom%edge_vect_n(i))*ut2_ad(i)
2223  ut2_ad(i) = 0.0_8
2224  ue_ad(i-1, npy, 1) = ue_ad(i-1, npy, 1) + geom%edge_vect_n(i)*ut1_ad(i)
2225  ue_ad(i, npy, 1) = ue_ad(i, npy, 1) + (1.-geom%edge_vect_n(i))*ut1_ad(i)
2226  ut1_ad(i) = 0.0_8
2227  end if
2228  end do
2229  end if
2230 
2231  if (js .eq. 1) then
2232  do i=ie,is,-1
2233  ut3_ad(i) = ut3_ad(i) + ue_ad(i, 1, 3)
2234  ue_ad(i, 1, 3) = 0.0_8
2235  ut2_ad(i) = ut2_ad(i) + ue_ad(i, 1, 2)
2236  ue_ad(i, 1, 2) = 0.0_8
2237  ut1_ad(i) = ut1_ad(i) + ue_ad(i, 1, 1)
2238  ue_ad(i, 1, 1) = 0.0_8
2239  end do
2240  do i=ie,is,-1
2241  if (i .le. im2) then
2242  ue_ad(i+1, 1, 3) = ue_ad(i+1, 1, 3) + geom%edge_vect_s(i)*ut3_ad(i)
2243  ue_ad(i, 1, 3) = ue_ad(i, 1, 3) + (1.-geom%edge_vect_s(i))*ut3_ad(i)
2244  ut3_ad(i) = 0.0_8
2245  ue_ad(i+1, 1, 2) = ue_ad(i+1, 1, 2) + geom%edge_vect_s(i)*ut2_ad(i)
2246  ue_ad(i, 1, 2) = ue_ad(i, 1, 2) + (1.-geom%edge_vect_s(i))*ut2_ad(i)
2247  ut2_ad(i) = 0.0_8
2248  ue_ad(i+1, 1, 1) = ue_ad(i+1, 1, 1) + geom%edge_vect_s(i)*ut1_ad(i)
2249  ue_ad(i, 1, 1) = ue_ad(i, 1, 1) + (1.-geom%edge_vect_s(i))*ut1_ad(i)
2250  ut1_ad(i) = 0.0_8
2251  else
2252  ue_ad(i-1, 1, 3) = ue_ad(i-1, 1, 3) + geom%edge_vect_s(i)*ut3_ad(i)
2253  ue_ad(i, 1, 3) = ue_ad(i, 1, 3) + (1.-geom%edge_vect_s(i))*ut3_ad(i)
2254  ut3_ad(i) = 0.0_8
2255  ue_ad(i-1, 1, 2) = ue_ad(i-1, 1, 2) + geom%edge_vect_s(i)*ut2_ad(i)
2256  ue_ad(i, 1, 2) = ue_ad(i, 1, 2) + (1.-geom%edge_vect_s(i))*ut2_ad(i)
2257  ut2_ad(i) = 0.0_8
2258  ue_ad(i-1, 1, 1) = ue_ad(i-1, 1, 1) + geom%edge_vect_s(i)*ut1_ad(i)
2259  ue_ad(i, 1, 1) = ue_ad(i, 1, 1) + (1.-geom%edge_vect_s(i))*ut1_ad(i)
2260  ut1_ad(i) = 0.0_8
2261  end if
2262  end do
2263  end if
2264 
2265  if (ie + 1 .eq. npx) then
2266  do j=je,js,-1
2267  vt3_ad(j) = vt3_ad(j) + ve_ad(npx, j, 3)
2268  ve_ad(npx, j, 3) = 0.0_8
2269  vt2_ad(j) = vt2_ad(j) + ve_ad(npx, j, 2)
2270  ve_ad(npx, j, 2) = 0.0_8
2271  vt1_ad(j) = vt1_ad(j) + ve_ad(npx, j, 1)
2272  ve_ad(npx, j, 1) = 0.0_8
2273  end do
2274  do j=je,js,-1
2275  if (j .le. jm2) then
2276  ve_ad(npx, j+1, 3) = ve_ad(npx, j+1, 3) + geom%edge_vect_e(j)*vt3_ad(j)
2277  ve_ad(npx, j, 3) = ve_ad(npx, j, 3) + (1.-geom%edge_vect_e(j))*vt3_ad(j)
2278  vt3_ad(j) = 0.0_8
2279  ve_ad(npx, j+1, 2) = ve_ad(npx, j+1, 2) + geom%edge_vect_e(j)*vt2_ad(j)
2280  ve_ad(npx, j, 2) = ve_ad(npx, j, 2) + (1.-geom%edge_vect_e(j))*vt2_ad(j)
2281  vt2_ad(j) = 0.0_8
2282  ve_ad(npx, j+1, 1) = ve_ad(npx, j+1, 1) + geom%edge_vect_e(j)*vt1_ad(j)
2283  ve_ad(npx, j, 1) = ve_ad(npx, j, 1) + (1.-geom%edge_vect_e(j))*vt1_ad(j)
2284  vt1_ad(j) = 0.0_8
2285  else
2286  ve_ad(npx, j-1, 3) = ve_ad(npx, j-1, 3) + geom%edge_vect_e(j)*vt3_ad(j)
2287  ve_ad(npx, j, 3) = ve_ad(npx, j, 3) + (1.-geom%edge_vect_e(j))*vt3_ad(j)
2288  vt3_ad(j) = 0.0_8
2289  ve_ad(npx, j-1, 2) = ve_ad(npx, j-1, 2) + geom%edge_vect_e(j)*vt2_ad(j)
2290  ve_ad(npx, j, 2) = ve_ad(npx, j, 2) + (1.-geom%edge_vect_e(j))*vt2_ad(j)
2291  vt2_ad(j) = 0.0_8
2292  ve_ad(npx, j-1, 1) = ve_ad(npx, j-1, 1) + geom%edge_vect_e(j)*vt1_ad(j)
2293  ve_ad(npx, j, 1) = ve_ad(npx, j, 1) + (1.-geom%edge_vect_e(j))*vt1_ad(j)
2294  vt1_ad(j) = 0.0_8
2295  end if
2296  end do
2297  end if
2298 
2299  if (is .eq. 1) then
2300  do j=je,js,-1
2301  vt3_ad(j) = vt3_ad(j) + ve_ad(1, j, 3)
2302  ve_ad(1, j, 3) = 0.0_8
2303  vt2_ad(j) = vt2_ad(j) + ve_ad(1, j, 2)
2304  ve_ad(1, j, 2) = 0.0_8
2305  vt1_ad(j) = vt1_ad(j) + ve_ad(1, j, 1)
2306  ve_ad(1, j, 1) = 0.0_8
2307  end do
2308  do j=je,js,-1
2309  if (j .le. jm2) then
2310  ve_ad(1, j+1, 3) = ve_ad(1, j+1, 3) + geom%edge_vect_w(j)*vt3_ad(j)
2311  ve_ad(1, j, 3) = ve_ad(1, j, 3) + (1.-geom%edge_vect_w(j))*vt3_ad(j)
2312  vt3_ad(j) = 0.0_8
2313  ve_ad(1, j+1, 2) = ve_ad(1, j+1, 2) + geom%edge_vect_w(j)*vt2_ad(j)
2314  ve_ad(1, j, 2) = ve_ad(1, j, 2) + (1.-geom%edge_vect_w(j))*vt2_ad(j)
2315  vt2_ad(j) = 0.0_8
2316  ve_ad(1, j+1, 1) = ve_ad(1, j+1, 1) + geom%edge_vect_w(j)*vt1_ad(j)
2317  ve_ad(1, j, 1) = ve_ad(1, j, 1) + (1.-geom%edge_vect_w(j))*vt1_ad(j)
2318  vt1_ad(j) = 0.0_8
2319  else
2320  ve_ad(1, j-1, 3) = ve_ad(1, j-1, 3) + geom%edge_vect_w(j)*vt3_ad(j)
2321  ve_ad(1, j, 3) = ve_ad(1, j, 3) + (1.-geom%edge_vect_w(j))*vt3_ad(j)
2322  vt3_ad(j) = 0.0_8
2323  ve_ad(1, j-1, 2) = ve_ad(1, j-1, 2) + geom%edge_vect_w(j)*vt2_ad(j)
2324  ve_ad(1, j, 2) = ve_ad(1, j, 2) + (1.-geom%edge_vect_w(j))*vt2_ad(j)
2325  vt2_ad(j) = 0.0_8
2326  ve_ad(1, j-1, 1) = ve_ad(1, j-1, 1) + geom%edge_vect_w(j)*vt1_ad(j)
2327  ve_ad(1, j, 1) = ve_ad(1, j, 1) + (1.-geom%edge_vect_w(j))*vt1_ad(j)
2328  vt1_ad(j) = 0.0_8
2329  end if
2330  end do
2331  end if
2332 
2333  do j=je+1,js-1,-1
2334  do i=ie+1,is,-1
2335  v3_ad(i-1, j, 3) = v3_ad(i-1, j, 3) + ve_ad(i, j, 3)
2336  v3_ad(i, j, 3) = v3_ad(i, j, 3) + ve_ad(i, j, 3)
2337  ve_ad(i, j, 3) = 0.0_8
2338  v3_ad(i-1, j, 2) = v3_ad(i-1, j, 2) + ve_ad(i, j, 2)
2339  v3_ad(i, j, 2) = v3_ad(i, j, 2) + ve_ad(i, j, 2)
2340  ve_ad(i, j, 2) = 0.0_8
2341  v3_ad(i-1, j, 1) = v3_ad(i-1, j, 1) + ve_ad(i, j, 1)
2342  v3_ad(i, j, 1) = v3_ad(i, j, 1) + ve_ad(i, j, 1)
2343  ve_ad(i, j, 1) = 0.0_8
2344  end do
2345  end do
2346 
2347  do j=je+1,js,-1
2348  do i=ie+1,is-1,-1
2349  v3_ad(i, j-1, 3) = v3_ad(i, j-1, 3) + ue_ad(i, j, 3)
2350  v3_ad(i, j, 3) = v3_ad(i, j, 3) + ue_ad(i, j, 3)
2351  ue_ad(i, j, 3) = 0.0_8
2352  v3_ad(i, j-1, 2) = v3_ad(i, j-1, 2) + ue_ad(i, j, 2)
2353  v3_ad(i, j, 2) = v3_ad(i, j, 2) + ue_ad(i, j, 2)
2354  ue_ad(i, j, 2) = 0.0_8
2355  v3_ad(i, j-1, 1) = v3_ad(i, j-1, 1) + ue_ad(i, j, 1)
2356  v3_ad(i, j, 1) = v3_ad(i, j, 1) + ue_ad(i, j, 1)
2357  ue_ad(i, j, 1) = 0.0_8
2358  end do
2359  end do
2360 
2361  do j=je+1,js-1,-1
2362  do i=ie+1,is-1,-1
2363  uatemp_ad(i, j, k) = uatemp_ad(i, j, k) + geom%vlon(i, j, 3)*v3_ad(i, j, 3)
2364  vatemp_ad(i, j, k) = vatemp_ad(i, j, k) + geom%vlat(i, j, 3)*v3_ad(i, j, 3)
2365  v3_ad(i, j, 3) = 0.0_8
2366  uatemp_ad(i, j, k) = uatemp_ad(i, j, k) + geom%vlon(i, j, 2)*v3_ad(i, j, 2)
2367  vatemp_ad(i, j, k) = vatemp_ad(i, j, k) + geom%vlat(i, j, 2)*v3_ad(i, j, 2)
2368  v3_ad(i, j, 2) = 0.0_8
2369  uatemp_ad(i, j, k) = uatemp_ad(i, j, k) + geom%vlon(i, j, 1)*v3_ad(i, j, 1)
2370  vatemp_ad(i, j, k) = vatemp_ad(i, j, k) + geom%vlat(i, j, 1)*v3_ad(i, j, 1)
2371  v3_ad(i, j, 1) = 0.0_8
2372  end do
2373  end do
2374  end do
2375 
2376  call mpp_update_domains_adm(vatemp, vatemp_ad, geom%domain, complete=.true.)
2377  call mpp_update_domains_adm(uatemp, uatemp_ad, geom%domain, complete=.true.)
2378 
2379  va_ad = va_ad + vatemp_ad(is:ie, js:je, :)
2380  ua_ad = ua_ad + uatemp_ad(is:ie, js:je, :)
2381 
2382 end subroutine a2d_ad
2383 
2384 ! ------------------------------------------------------------------------------
2385 
2386  subroutine d2a(geom, u, v, ua, va)
2387 
2388  !c2l_ord4
2389 
2390  !Called using model fields so defined on isd:ied,jsd:jed
2391 
2392  use mpp_domains_mod, only: mpp_update_domains, dgrid_ne
2393 
2394  implicit none
2395 
2396  type(fv3jedi_geom), intent(inout) :: geom
2397 
2398  real(kind=kind_real), intent(inout):: u(geom%isd:geom%ied ,geom%jsd:geom%jed+1,geom%npz)
2399  real(kind=kind_real), intent(inout):: v(geom%isd:geom%ied+1,geom%jsd:geom%jed ,geom%npz)
2400  real(kind=kind_real), intent(inout):: ua(geom%isd:geom%ied ,geom%jsd:geom%jed ,geom%npz)
2401  real(kind=kind_real), intent(inout):: va(geom%isd:geom%ied ,geom%jsd:geom%jed ,geom%npz)
2402 
2403 ! Local
2404 ! 4-pt Lagrange interpolation
2405  real(kind=kind_real) :: c1 = 1.125
2406  real(kind=kind_real) :: c2 = -0.125
2407  real(kind=kind_real) :: utmp(geom%isc:geom%iec, geom%jsc:geom%jec+1)
2408  real(kind=kind_real) :: vtmp(geom%isc:geom%iec+1,geom%jsc:geom%jec)
2409  real(kind=kind_real) :: wu(geom%isc:geom%iec, geom%jsc:geom%jec+1)
2410  real(kind=kind_real) :: wv(geom%isc:geom%iec+1,geom%jsc:geom%jec)
2411 
2412  integer i, j, k
2413  integer :: is, ie, js, je, npx, npy, npz
2414 
2415  is = geom%isc
2416  ie = geom%iec
2417  js = geom%jsc
2418  je = geom%jec
2419  npx = geom%npx
2420  npy = geom%npy
2421  npz = geom%npz
2422 
2423  call mpp_update_domains(u, v, geom%domain, gridtype=dgrid_ne)
2424 
2425 !$OMP parallel do default(none) shared(is,ie,js,je,npz,npx,npy,c2,c1, &
2426 !$OMP u,v,ua,va) &
2427 !$OMP private(utmp, vtmp, wu, wv)
2428  do k=1,npz
2429 
2430  do j=max(2,js),min(npy-2,je)
2431  do i=max(2,is),min(npx-2,ie)
2432  utmp(i,j) = c2*(u(i,j-1,k)+u(i,j+2,k)) + c1*(u(i,j,k)+u(i,j+1,k))
2433  vtmp(i,j) = c2*(v(i-1,j,k)+v(i+2,j,k)) + c1*(v(i,j,k)+v(i+1,j,k))
2434  enddo
2435  enddo
2436 
2437  if ( js==1 ) then
2438  do i=is,ie+1
2439  wv(i,1) = v(i,1,k)*geom%dy(i,1)
2440  enddo
2441  do i=is,ie
2442  vtmp(i,1) = 2.*(wv(i,1) + wv(i+1,1)) / (geom%dy(i,1)+geom%dy(i+1,1))
2443  utmp(i,1) = 2.*(u(i,1,k)*geom%dx(i,1) + u(i,2,k)*geom%dx(i,2)) &
2444  / ( geom%dx(i,1) + geom%dx(i,2))
2445  enddo
2446  endif
2447 
2448  if ( (je+1)==npy ) then
2449  j = npy-1
2450  do i=is,ie+1
2451  wv(i,j) = v(i,j,k)*geom%dy(i,j)
2452  enddo
2453  do i=is,ie
2454  vtmp(i,j) = 2.*(wv(i,j) + wv(i+1,j)) / (geom%dy(i,j)+geom%dy(i+1,j))
2455  utmp(i,j) = 2.*(u(i,j,k)*geom%dx(i,j) + u(i,j+1,k)*geom%dx(i,j+1)) &
2456  / ( geom%dx(i,j) + geom%dx(i,j+1))
2457  enddo
2458  endif
2459 
2460  if ( is==1 ) then
2461  i = 1
2462  do j=js,je
2463  wv(1,j) = v(1,j,k)*geom%dy(1,j)
2464  wv(2,j) = v(2,j,k)*geom%dy(2,j)
2465  enddo
2466  do j=js,je+1
2467  wu(i,j) = u(i,j,k)*geom%dx(i,j)
2468  enddo
2469  do j=js,je
2470  utmp(i,j) = 2.*(wu(i,j) + wu(i,j+1))/(geom%dx(i,j)+geom%dx(i,j+1))
2471  vtmp(i,j) = 2.*(wv(1,j) + wv(2,j ))/(geom%dy(1,j)+geom%dy(2,j))
2472  enddo
2473  endif
2474 
2475  if ( (ie+1)==npx) then
2476  i = npx-1
2477  do j=js,je
2478  wv(i, j) = v(i, j,k)*geom%dy(i, j)
2479  wv(i+1,j) = v(i+1,j,k)*geom%dy(i+1,j)
2480  enddo
2481  do j=js,je+1
2482  wu(i,j) = u(i,j,k)*geom%dx(i,j)
2483  enddo
2484  do j=js,je
2485  utmp(i,j) = 2.*(wu(i,j) + wu(i, j+1))/(geom%dx(i,j)+geom%dx(i,j+1))
2486  vtmp(i,j) = 2.*(wv(i,j) + wv(i+1,j ))/(geom%dy(i,j)+geom%dy(i+1,j))
2487  enddo
2488  endif
2489 
2490  !Transform local a-grid winds into latitude-longitude coordinates
2491  do j=js,je
2492  do i=is,ie
2493  ua(i,j,k) = geom%a11(i,j)*utmp(i,j) + geom%a12(i,j)*vtmp(i,j)
2494  va(i,j,k) = geom%a21(i,j)*utmp(i,j) + geom%a22(i,j)*vtmp(i,j)
2495  enddo
2496  enddo
2497 
2498  enddo
2499 
2500 end subroutine d2a
2501 
2502 ! ------------------------------------------------------------------------------
2503 
2504  SUBROUTINE d2a_ad(geom, u_ad, v_ad, ua_ad, va_ad)
2505 
2506  use fv_mp_adm_mod, only: mpp_update_domains_adm
2507  use mpp_domains_mod, only: dgrid_ne
2508 
2509  IMPLICIT NONE
2510  TYPE(fv3jedi_geom), INTENT(INOUT) :: geom
2511  REAL(kind=kind_real) :: u(geom%isd:geom%ied, geom%jsd&
2512 & :geom%jed+1, geom%npz)
2513  REAL(kind=kind_real), INTENT(INOUT) :: u_ad(geom%isd:geom%ied, geom%&
2514 & jsd:geom%jed+1, geom%npz)
2515  REAL(kind=kind_real) :: v(geom%isd:geom%ied+1, geom%&
2516 & jsd:geom%jed, geom%npz)
2517  REAL(kind=kind_real), INTENT(INOUT) :: v_ad(geom%isd:geom%ied+1, &
2518 & geom%jsd:geom%jed, geom%npz)
2519  REAL(kind=kind_real) :: ua(geom%isd:geom%ied, geom%&
2520 & jsd:geom%jed, geom%npz)
2521  REAL(kind=kind_real), INTENT(INOUT) :: ua_ad(geom%isd:geom%ied, geom&
2522 & %jsd:geom%jed, geom%npz)
2523  REAL(kind=kind_real) :: va(geom%isd:geom%ied, geom%&
2524 & jsd:geom%jed, geom%npz)
2525  REAL(kind=kind_real), INTENT(INOUT) :: va_ad(geom%isd:geom%ied, geom&
2526 & %jsd:geom%jed, geom%npz)
2527 ! Local
2528 ! 4-pt Lagrange interpolation
2529  REAL(kind=kind_real), SAVE :: c1=1.125
2530  REAL(kind=kind_real), SAVE :: c2=-0.125
2531  REAL(kind=kind_real) :: utmp(geom%isc:geom%iec, geom%jsc:geom%jec+1)
2532  REAL(kind=kind_real) :: utmp_ad(geom%isc:geom%iec, geom%jsc:geom%jec&
2533 & +1)
2534  REAL(kind=kind_real) :: vtmp(geom%isc:geom%iec+1, geom%jsc:geom%jec)
2535  REAL(kind=kind_real) :: vtmp_ad(geom%isc:geom%iec+1, geom%jsc:geom%&
2536 & jec)
2537  REAL(kind=kind_real) :: wu(geom%isc:geom%iec, geom%jsc:geom%jec+1)
2538  REAL(kind=kind_real) :: wu_ad(geom%isc:geom%iec, geom%jsc:geom%jec+1&
2539 & )
2540  REAL(kind=kind_real) :: wv(geom%isc:geom%iec+1, geom%jsc:geom%jec)
2541  REAL(kind=kind_real) :: wv_ad(geom%isc:geom%iec+1, geom%jsc:geom%jec&
2542 & )
2543  INTEGER :: i, j, k
2544  INTEGER :: is, ie, js, je, npx, npy, npz
2545  INTRINSIC max
2546  INTRINSIC min
2547  INTEGER :: max1
2548  INTEGER :: max2
2549  INTEGER :: min1
2550  INTEGER :: min2
2551  REAL(kind=kind_real) :: temp_ad
2552  REAL(kind=kind_real) :: temp_ad0
2553  REAL(kind=kind_real) :: temp_ad1
2554  REAL(kind=kind_real) :: temp_ad2
2555  REAL(kind=kind_real) :: temp_ad3
2556  REAL(kind=kind_real) :: temp_ad4
2557  REAL(kind=kind_real) :: temp_ad5
2558  REAL(kind=kind_real) :: temp_ad6
2559  INTEGER :: ad_from
2560  INTEGER :: ad_to
2561  INTEGER :: ad_from0
2562  INTEGER :: ad_to0
2563  INTEGER :: branch
2564  is = geom%isc
2565  ie = geom%iec
2566  js = geom%jsc
2567  je = geom%jec
2568  npx = geom%npx
2569  npy = geom%npy
2570  npz = geom%npz
2571 !$OMP parallel do default(none) shared(is,ie,js,je,npz,npx,npy,c2,c1, &
2572 !$OMP u,v,ua,va) &
2573 !$OMP private(utmp, vtmp, wu, wv)
2574  DO k=1,npz
2575  IF (2 .LT. js) THEN
2576  max1 = js
2577  ELSE
2578  max1 = 2
2579  END IF
2580  IF (npy - 2 .GT. je) THEN
2581  min1 = je
2582  ELSE
2583  min1 = npy - 2
2584  END IF
2585  ad_from0 = max1
2586  DO j=ad_from0,min1
2587  IF (2 .LT. is) THEN
2588  max2 = is
2589  ELSE
2590  max2 = 2
2591  END IF
2592  IF (npx - 2 .GT. ie) THEN
2593  min2 = ie
2594  ELSE
2595  min2 = npx - 2
2596  END IF
2597  ad_from = max2
2598  CALL pushinteger4(i)
2599  i = min2 + 1
2600  CALL pushinteger4(i - 1)
2601  CALL pushinteger4(ad_from)
2602  END DO
2603  CALL pushinteger4(j - 1)
2604  CALL pushinteger4(ad_from0)
2605  IF (js .EQ. 1) THEN
2606  CALL pushinteger4(i)
2607  CALL pushcontrol1b(0)
2608  ELSE
2609  CALL pushcontrol1b(1)
2610  END IF
2611  IF (je + 1 .EQ. npy) THEN
2612  j = npy - 1
2613  CALL pushinteger4(i)
2614  CALL pushcontrol1b(0)
2615  ELSE
2616  CALL pushcontrol1b(1)
2617  END IF
2618  IF (is .EQ. 1) THEN
2619  CALL pushinteger4(i)
2620  i = 1
2621  CALL pushinteger4(j)
2622  CALL pushcontrol1b(0)
2623  ELSE
2624  CALL pushcontrol1b(1)
2625  END IF
2626  IF (ie + 1 .EQ. npx) THEN
2627  CALL pushinteger4(i)
2628  i = npx - 1
2629  CALL pushinteger4(j)
2630  CALL pushcontrol1b(1)
2631  ELSE
2632  CALL pushcontrol1b(0)
2633  END IF
2634  CALL pushinteger4(j)
2635 !Transform local a-grid winds into latitude-longitude coordinates
2636  DO j=js,je
2637  CALL pushinteger4(i)
2638  END DO
2639  END DO
2640  vtmp_ad = 0.0_8
2641  wu_ad = 0.0_8
2642  wv_ad = 0.0_8
2643  utmp_ad = 0.0_8
2644  DO k=npz,1,-1
2645  DO j=je,js,-1
2646  DO i=ie,is,-1
2647  utmp_ad(i, j) = utmp_ad(i, j) + geom%a11(i, j)*ua_ad(i, j, k) &
2648 & + geom%a21(i, j)*va_ad(i, j, k)
2649  vtmp_ad(i, j) = vtmp_ad(i, j) + geom%a12(i, j)*ua_ad(i, j, k) &
2650 & + geom%a22(i, j)*va_ad(i, j, k)
2651  va_ad(i, j, k) = 0.0_8
2652  ua_ad(i, j, k) = 0.0_8
2653  END DO
2654  CALL popinteger4(i)
2655  END DO
2656  CALL popinteger4(j)
2657  CALL popcontrol1b(branch)
2658  IF (branch .NE. 0) THEN
2659  DO j=je,js,-1
2660  temp_ad5 = 2.*vtmp_ad(i, j)/(geom%dy(i, j)+geom%dy(i+1, j))
2661  wv_ad(i, j) = wv_ad(i, j) + temp_ad5
2662  wv_ad(i+1, j) = wv_ad(i+1, j) + temp_ad5
2663  vtmp_ad(i, j) = 0.0_8
2664  temp_ad6 = 2.*utmp_ad(i, j)/(geom%dx(i, j)+geom%dx(i, j+1))
2665  wu_ad(i, j) = wu_ad(i, j) + temp_ad6
2666  wu_ad(i, j+1) = wu_ad(i, j+1) + temp_ad6
2667  utmp_ad(i, j) = 0.0_8
2668  END DO
2669  DO j=je+1,js,-1
2670  u_ad(i, j, k) = u_ad(i, j, k) + geom%dx(i, j)*wu_ad(i, j)
2671  wu_ad(i, j) = 0.0_8
2672  END DO
2673  DO j=je,js,-1
2674  v_ad(i+1, j, k) = v_ad(i+1, j, k) + geom%dy(i+1, j)*wv_ad(i+1&
2675 & , j)
2676  wv_ad(i+1, j) = 0.0_8
2677  v_ad(i, j, k) = v_ad(i, j, k) + geom%dy(i, j)*wv_ad(i, j)
2678  wv_ad(i, j) = 0.0_8
2679  END DO
2680  CALL popinteger4(j)
2681  CALL popinteger4(i)
2682  END IF
2683  CALL popcontrol1b(branch)
2684  IF (branch .EQ. 0) THEN
2685  DO j=je,js,-1
2686  temp_ad3 = 2.*vtmp_ad(i, j)/(geom%dy(1, j)+geom%dy(2, j))
2687  wv_ad(1, j) = wv_ad(1, j) + temp_ad3
2688  wv_ad(2, j) = wv_ad(2, j) + temp_ad3
2689  vtmp_ad(i, j) = 0.0_8
2690  temp_ad4 = 2.*utmp_ad(i, j)/(geom%dx(i, j)+geom%dx(i, j+1))
2691  wu_ad(i, j) = wu_ad(i, j) + temp_ad4
2692  wu_ad(i, j+1) = wu_ad(i, j+1) + temp_ad4
2693  utmp_ad(i, j) = 0.0_8
2694  END DO
2695  DO j=je+1,js,-1
2696  u_ad(i, j, k) = u_ad(i, j, k) + geom%dx(i, j)*wu_ad(i, j)
2697  wu_ad(i, j) = 0.0_8
2698  END DO
2699  DO j=je,js,-1
2700  v_ad(2, j, k) = v_ad(2, j, k) + geom%dy(2, j)*wv_ad(2, j)
2701  wv_ad(2, j) = 0.0_8
2702  v_ad(1, j, k) = v_ad(1, j, k) + geom%dy(1, j)*wv_ad(1, j)
2703  wv_ad(1, j) = 0.0_8
2704  END DO
2705  CALL popinteger4(j)
2706  CALL popinteger4(i)
2707  END IF
2708  CALL popcontrol1b(branch)
2709  IF (branch .EQ. 0) THEN
2710  DO i=ie,is,-1
2711  temp_ad1 = 2.*utmp_ad(i, j)/(geom%dx(i, j)+geom%dx(i, j+1))
2712  u_ad(i, j, k) = u_ad(i, j, k) + geom%dx(i, j)*temp_ad1
2713  u_ad(i, j+1, k) = u_ad(i, j+1, k) + geom%dx(i, j+1)*temp_ad1
2714  utmp_ad(i, j) = 0.0_8
2715  temp_ad2 = 2.*vtmp_ad(i, j)/(geom%dy(i, j)+geom%dy(i+1, j))
2716  wv_ad(i, j) = wv_ad(i, j) + temp_ad2
2717  wv_ad(i+1, j) = wv_ad(i+1, j) + temp_ad2
2718  vtmp_ad(i, j) = 0.0_8
2719  END DO
2720  DO i=ie+1,is,-1
2721  v_ad(i, j, k) = v_ad(i, j, k) + geom%dy(i, j)*wv_ad(i, j)
2722  wv_ad(i, j) = 0.0_8
2723  END DO
2724  CALL popinteger4(i)
2725  END IF
2726  CALL popcontrol1b(branch)
2727  IF (branch .EQ. 0) THEN
2728  DO i=ie,is,-1
2729  temp_ad = 2.*utmp_ad(i, 1)/(geom%dx(i, 1)+geom%dx(i, 2))
2730  u_ad(i, 1, k) = u_ad(i, 1, k) + geom%dx(i, 1)*temp_ad
2731  u_ad(i, 2, k) = u_ad(i, 2, k) + geom%dx(i, 2)*temp_ad
2732  utmp_ad(i, 1) = 0.0_8
2733  temp_ad0 = 2.*vtmp_ad(i, 1)/(geom%dy(i, 1)+geom%dy(i+1, 1))
2734  wv_ad(i, 1) = wv_ad(i, 1) + temp_ad0
2735  wv_ad(i+1, 1) = wv_ad(i+1, 1) + temp_ad0
2736  vtmp_ad(i, 1) = 0.0_8
2737  END DO
2738  DO i=ie+1,is,-1
2739  v_ad(i, 1, k) = v_ad(i, 1, k) + geom%dy(i, 1)*wv_ad(i, 1)
2740  wv_ad(i, 1) = 0.0_8
2741  END DO
2742  CALL popinteger4(i)
2743  END IF
2744  CALL popinteger4(ad_from0)
2745  CALL popinteger4(ad_to0)
2746  DO j=ad_to0,ad_from0,-1
2747  CALL popinteger4(ad_from)
2748  CALL popinteger4(ad_to)
2749  DO i=ad_to,ad_from,-1
2750  v_ad(i-1, j, k) = v_ad(i-1, j, k) + c2*vtmp_ad(i, j)
2751  v_ad(i+2, j, k) = v_ad(i+2, j, k) + c2*vtmp_ad(i, j)
2752  v_ad(i, j, k) = v_ad(i, j, k) + c1*vtmp_ad(i, j)
2753  v_ad(i+1, j, k) = v_ad(i+1, j, k) + c1*vtmp_ad(i, j)
2754  vtmp_ad(i, j) = 0.0_8
2755  u_ad(i, j-1, k) = u_ad(i, j-1, k) + c2*utmp_ad(i, j)
2756  u_ad(i, j+2, k) = u_ad(i, j+2, k) + c2*utmp_ad(i, j)
2757  u_ad(i, j, k) = u_ad(i, j, k) + c1*utmp_ad(i, j)
2758  u_ad(i, j+1, k) = u_ad(i, j+1, k) + c1*utmp_ad(i, j)
2759  utmp_ad(i, j) = 0.0_8
2760  END DO
2761  CALL popinteger4(i)
2762  END DO
2763  END DO
2764  CALL mpp_update_domains_adm(u, u_ad, v, v_ad, geom%domain, gridtype=&
2765 & dgrid_ne)
2766  END SUBROUTINE d2a_ad
2767 
2768 ! ------------------------------------------------------------------------------
2769 
2770 end module wind_vt_mod
fv_mp_adm_mod::mpp_update_domains_adm
Definition: fv_mp_adm.F90:91
wind_vt_mod::extrap_corner_adm
subroutine, public extrap_corner_adm(p0, p1, p2, q1, q1_ad, q2, q2_ad, extrap_corner_ad)
Definition: wind_variables_mod.f90:1932
popinteger4
subroutine popinteger4(x)
Definition: adBuffer.f:541
wind_vt_mod::psichi_to_udvd_adm
subroutine, public psichi_to_udvd_adm(geom, psi_ad, chi_ad, u_ad, v_ad)
Definition: wind_variables_mod.f90:1152
wind_vt_mod::psichi_to_uava_adm
subroutine, public psichi_to_uava_adm(geom, psi_ad, chi_ad, ua_ad, va_ad)
Definition: wind_variables_mod.f90:345
fv3jedi_constants_mod::rad2deg
real(kind=kind_real), parameter, public rad2deg
Definition: fv3jedi_constants_mod.f90:13
mpp_domains_mod::mpp_update_domains
Definition: mpp_domains.F90:1068
wind_vt_mod::vortdivg_to_psichi
subroutine, public vortdivg_to_psichi(geom, vort, divg, psi, chi)
Definition: wind_variables_mod.f90:204
wind_vt_mod::great_circle_dist
real(kind=kind_real) function great_circle_dist(q1, q2)
Definition: wind_variables_mod.f90:1132
wind_vt_mod::a2b_ord4_adm
subroutine, public a2b_ord4_adm(qin, qin_ad, qout, qout_ad, geom, npx, npy, is, ie, js, je, ng)
Definition: wind_variables_mod.f90:1239
wind_vt_mod::a2b_ord4
subroutine, public a2b_ord4(qin, qout, geom, npx, npy, is, ie, js, je, ng)
Definition: wind_variables_mod.f90:926
wind_vt_mod::d2a
subroutine, public d2a(geom, u, v, ua, va)
Definition: wind_variables_mod.f90:2387
pushcontrol1b
subroutine pushcontrol1b(cc)
Definition: adBuffer.f:115
fv3jedi_geom_mod::fv3jedi_geom
Fortran derived type to hold geometry data for the FV3JEDI model.
Definition: fv3jedi_geom_mod.f90:34
wind_vt_mod::a2d_ad
subroutine, public a2d_ad(geom, ua_ad, va_ad, ud_ad, vd_ad)
Definition: wind_variables_mod.f90:2123
mpp_domains_mod
Definition: mpp_domains.F90:128
wind_vt_mod::sfc_10m_winds
subroutine, public sfc_10m_winds(geom, usrf, vsrf, f10r, spd10m, dir10m)
Definition: wind_variables_mod.f90:44
fv3jedi_constants_mod
Definition: fv3jedi_constants_mod.f90:6
wind_vt_mod::uv_to_vortdivg
subroutine, public uv_to_vortdivg(geom, u, v, ua, va, vort, divg)
Definition: wind_variables_mod.f90:105
wind_vt_mod::extrap_corner
real(kind=kind_real) function extrap_corner(p0, p1, p2, q1, q2)
Definition: wind_variables_mod.f90:1116
wind_vt_mod::a2d
subroutine, public a2d(geom, ua, va, ud, vd)
Definition: wind_variables_mod.f90:1950
fv_mp_adm_mod
Definition: fv_mp_adm.F90:1
wind_vt_mod::psichi_to_uava
subroutine, public psichi_to_uava(geom, psi, chi, ua, va)
Definition: wind_variables_mod.f90:301
mpp_domains_mod::mpp_get_boundary
Definition: mpp_domains.F90:1900
wind_vt_mod::gauss_seidel
subroutine, public gauss_seidel(geom, x, b)
Definition: wind_variables_mod.f90:239
wind_vt_mod
Variable transforms on wind variables for fv3-jedi Daniel Holdaway, NASA/JCSDA.
Definition: wind_variables_mod.f90:9
popcontrol1b
subroutine popcontrol1b(cc)
Definition: adBuffer.f:120
wind_vt_mod::edge_interpolate4
real(kind=kind_real) function edge_interpolate4(ua, dxa)
Definition: wind_variables_mod.f90:828
fv3jedi_geom_mod
Fortran module handling geometry for the FV3 model.
Definition: fv3jedi_geom_mod.f90:8
max
#define max(a, b)
Definition: mosaic_util.h:33
fv3jedi_kinds_mod
Definition: fv3jedi_kinds_mod.f90:6
wind_vt_mod::d2a_ad
subroutine, public d2a_ad(geom, u_ad, v_ad, ua_ad, va_ad)
Definition: wind_variables_mod.f90:2505
wind_vt_mod::d_to_ac
subroutine, public d_to_ac(geom, u, v, ua, va, uc, vc)
Definition: wind_variables_mod.f90:457
fv3jedi_kinds_mod::kind_real
integer, parameter, public kind_real
Definition: fv3jedi_kinds_mod.f90:13
min
#define min(a, b)
Definition: mosaic_util.h:32
wind_vt_mod::psichi_to_udvd
subroutine, public psichi_to_udvd(geom, psi, chi, u, v)
Definition: wind_variables_mod.f90:404
wind_vt_mod::d2a2c_vect
subroutine, public d2a2c_vect(geom, u, v, ua, va, uc, vc, ut, vt, dord4)
Definition: wind_variables_mod.f90:534
wind_vt_mod::divergence_corner
subroutine, public divergence_corner(geom, u, v, ua, va, divg_d)
Definition: wind_variables_mod.f90:846
pushinteger4
subroutine pushinteger4(x)
Definition: adBuffer.f:484
crtm_parameters::pi
real(fp), parameter, public pi
Definition: CRTM_Parameters.f90:86