Line data Source code
1 : !-----------------------------------------------------------------------
2 : !BOP
3 : ! !IROUTINE: p_d_adjust --- complete full physics update
4 : !
5 : ! !INTERFACE:
6 14592 : subroutine p_d_adjust(grid, tracer, pelnxy, pkxy, pkzxy, zvir, &
7 14592 : cap3v, delpxy, ptxy, pexy, psxy, ptop)
8 :
9 : ! !USES:
10 : use shr_kind_mod, only: r8 => shr_kind_r8
11 : use dynamics_vars, only : T_FVDYCORE_GRID
12 : #if defined( SPMD )
13 : use parutilitiesmodule, only: parcollective, sumop
14 : #endif
15 : use shr_reprosum_mod, only : shr_reprosum_calc, shr_reprosum_tolExceeded, &
16 : shr_reprosum_reldiffmax, &
17 : shr_reprosum_recompute
18 : use cam_logfile, only : iulog
19 : use perf_mod
20 :
21 : !-----------------------------------------------------------------------
22 : implicit none
23 :
24 : ! !INPUT PARAMETERS:
25 : type (T_FVDYCORE_GRID), intent(in) :: grid
26 : real(r8), intent(in) :: zvir
27 : real(r8), intent(in) :: ptop
28 : real(r8), intent(in) :: cap3v( grid%ifirstxy:grid%ilastxy, &
29 : grid%jfirstxy:grid%jlastxy, grid%km) ! cappa
30 :
31 : ! !INPUT/OUTPUT PARAMETERS:
32 : real(r8), intent(inout) :: tracer(grid%ifirstxy:grid%ilastxy, &
33 : grid%jfirstxy:grid%jlastxy,grid%km,grid%ntotq) ! constituents
34 : real(r8), intent(inout) :: delpxy(grid%ifirstxy:grid%ilastxy, &
35 : grid%jfirstxy:grid%jlastxy,grid%km) ! pressure difference
36 : real(r8), intent(inout) :: ptxy (grid%ifirstxy:grid%ilastxy, &
37 : grid%jfirstxy:grid%jlastxy, grid%km) ! Virtual pot temp
38 : real(r8), intent(inout) :: pexy(grid%ifirstxy:grid%ilastxy, &
39 : grid%km+1,grid%jfirstxy:grid%jlastxy)
40 :
41 : ! !OUTPUT PARAMETERS
42 : real(r8), intent(out) :: psxy(grid%ifirstxy:grid%ilastxy, grid%jfirstxy:grid%jlastxy) ! surf. press
43 : real(r8), intent(out) :: pelnxy(grid%ifirstxy:grid%ilastxy,grid%km+1,grid%jfirstxy:grid%jlastxy) ! interface pres
44 : real(r8), intent(out) :: pkzxy(grid%ifirstxy:grid%ilastxy,grid%jfirstxy:grid%jlastxy, grid%km) ! Layer-mean value of PK
45 : real(r8), intent(out) :: pkxy(grid%ifirstxy:grid%ilastxy,grid%jfirstxy:grid%jlastxy, grid%km+1) ! PE**cappa
46 :
47 : ! !DESCRIPTION:
48 : !
49 : ! Complete adjustment of quantities after physics update
50 : !
51 : ! !REVISION HISTORY:
52 : ! 00.06.01 Grant? Creation
53 : ! 01.06.08 AAM Created from p_d_coupling
54 : ! 02.04.24 WS New mod_comm interface
55 : ! 02.05.01 WS Fix of S.-J. and Phil to peln, pk update
56 : ! 03.03.31 BAB dry mass adjustment moved to dme_adjust, just finish up here
57 : ! 05.07.06 WS Use grid argument to get all grid-related data
58 : ! 05.09.23 WS Transitioned to XY variables only
59 : ! 06.07.01 WS Transitioned tracers q3 to T_TRACERS
60 : ! 08.06.25 PW Added call to fixed point reproducible sum
61 : !
62 : !EOP
63 : !-----------------------------------------------------------------------
64 : !BOC
65 : ! !LOCAL VARIABLES:
66 : real(r8), parameter :: D0_0 = 0.0_r8
67 : real(r8), parameter :: D1_0 = 1.0_r8
68 :
69 29184 : real(r8) :: pole(grid%ifirstxy:grid%ilastxy,grid%km,grid%ntotq+2)
70 : ! Array containing local pole values
71 29184 : real(r8) :: pole_sum(grid%km,grid%ntotq+2) ! Array containing average of all pole values
72 29184 : real(r8) :: rel_diff(2,grid%km,grid%ntotq+2)
73 :
74 : real(r8) :: &
75 29184 : cap3vi(grid%ifirstxy:grid%ilastxy,grid%jfirstxy:grid%jlastxy,grid%km+1)
76 : real(r8) :: &
77 29184 : pkln(grid%ifirstxy:grid%ilastxy,grid%km+1,grid%jfirstxy:grid%jlastxy)
78 :
79 14592 : real(r8),allocatable :: pole_tmp(:)
80 :
81 : integer :: i, k, m, j ! indices
82 : integer :: im, jm, km, ntotq, lim
83 : integer :: ifirstxy, ilastxy, jfirstxy, jlastxy
84 :
85 : logical :: write_warning
86 : logical :: high_alt
87 :
88 : !---------------------------End Local workspace-------------------------
89 :
90 : !
91 : ! ----------------------------------------------------
92 : ! Complete update of dynamics variables
93 : ! ----------------------------------------------------
94 : !
95 14592 : high_alt = grid%high_alt
96 14592 : im = grid%im
97 14592 : jm = grid%jm
98 14592 : km = grid%km
99 14592 : ntotq = grid%ntotq
100 :
101 14592 : ifirstxy = grid%ifirstxy
102 14592 : ilastxy = grid%ilastxy
103 14592 : jfirstxy = grid%jfirstxy
104 14592 : jlastxy = grid%jlastxy
105 :
106 14592 : lim = (ilastxy-ifirstxy) + 1
107 :
108 : ! Average the pole values (WS 2006/02/16, bug fix)
109 :
110 14592 : if (jfirstxy==1) then
111 :
112 : !$omp parallel do private(i,k,m)
113 7524 : do k = 1, km
114 182400 : do i = ifirstxy, ilastxy
115 175104 : pole(i,k,1) = delpxy(i,1,k)
116 182400 : pole(i,k,2) = ptxy(i,1,k)
117 : enddo
118 2867556 : do m = 1, ntotq
119 71508096 : do i = ifirstxy, ilastxy
120 71500800 : pole(i,k,m+2) = tracer(i,1,k,m)
121 : enddo
122 : enddo
123 : enddo
124 :
125 228 : call t_startf("pdadj_reprosum")
126 : call shr_reprosum_calc(pole, pole_sum, lim, lim, km*(ntotq+2), gbl_count=im, &
127 228 : commid=grid%commxy_x, rel_diff=rel_diff) ! South pole
128 228 : call t_stopf("pdadj_reprosum")
129 :
130 : ! check that "fast" reproducible sum is accurate enough. If not, calculate
131 : ! using old method
132 228 : write_warning = .false.
133 228 : if (grid%myidxy_x == 0) write_warning = .true.
134 228 : if ( shr_reprosum_tolExceeded('p_d_adjust/South Pole', km*(ntotq+2), &
135 : write_warning, iulog, rel_diff) ) then
136 0 : if ( shr_reprosum_recompute ) then
137 0 : call t_startf("pdadj_sumfix")
138 0 : allocate( pole_tmp(im) )
139 0 : do m = 1, ntotq+2
140 0 : do k = 1, km
141 0 : if (rel_diff(1,k,m) > shr_reprosum_reldiffmax) then
142 0 : pole_tmp(:) = D0_0
143 0 : do i = ifirstxy, ilastxy
144 0 : pole_tmp(i) = pole(i,k,m)
145 : enddo
146 : #if defined(SPMD)
147 0 : call parcollective(grid%commxy_x,sumop,im,pole_tmp)
148 : #endif
149 0 : pole_sum(k,m) = D0_0
150 0 : do i = 1, im
151 0 : pole_sum(k,m) = pole_sum(k,m) + pole_tmp(i)
152 : enddo
153 : endif
154 : enddo
155 : enddo
156 0 : deallocate( pole_tmp )
157 0 : call t_stopf("pdadj_sumfix")
158 : endif
159 : endif
160 :
161 : ! save results
162 : !$omp parallel do private(i,k,m)
163 7524 : do k = 1, km
164 : ! normalize first
165 2881920 : do m = 1, ntotq+2
166 2881920 : pole_sum(k,m) = pole_sum(k,m)/im
167 : enddo
168 182400 : do i = ifirstxy,ilastxy
169 175104 : delpxy(i,1,k) = pole_sum(k,1)
170 182400 : ptxy(i,1,k) = pole_sum(k,2)
171 : enddo
172 2867556 : do m = 1, ntotq
173 71508096 : do i = ifirstxy,ilastxy
174 71500800 : tracer(i,1,k,m) = pole_sum(k,m+2)
175 : enddo
176 : enddo
177 : enddo
178 :
179 : endif ! jfirstxy==1
180 :
181 14592 : if (jlastxy==jm) then
182 :
183 : !$omp parallel do private(i,k,m)
184 7524 : do k = 1, km
185 182400 : do i = ifirstxy, ilastxy
186 175104 : pole(i,k,1) = delpxy(i,jm,k)
187 182400 : pole(i,k,2) = ptxy(i,jm,k)
188 : enddo
189 2867556 : do m = 1, ntotq
190 71508096 : do i = ifirstxy, ilastxy
191 71500800 : pole(i,k,m+2) = tracer(i,jm,k,m)
192 : enddo
193 : enddo
194 : enddo
195 :
196 228 : call t_startf("pdadj_reprosum")
197 : call shr_reprosum_calc(pole, pole_sum, lim, lim, km*(ntotq+2), gbl_count=im, &
198 228 : commid=grid%commxy_x, rel_diff=rel_diff) ! North pole
199 228 : call t_stopf("pdadj_reprosum")
200 :
201 : ! check that "fast" reproducible sum is accurate enough. If not, calculate
202 : ! using old method
203 228 : write_warning = .false.
204 228 : if (grid%myidxy_x == 0) write_warning = .true.
205 228 : if ( shr_reprosum_tolExceeded('p_d_adjust/Nouth Pole', km*(ntotq+2), &
206 : write_warning, iulog, rel_diff) ) then
207 0 : if ( shr_reprosum_recompute ) then
208 0 : call t_startf("pdadj_sumfix")
209 0 : allocate( pole_tmp(im) )
210 0 : do m = 1, ntotq+2
211 0 : do k = 1, km
212 0 : if (rel_diff(1,k,m) > shr_reprosum_reldiffmax) then
213 0 : pole_tmp(:) = D0_0
214 0 : do i = ifirstxy, ilastxy
215 0 : pole_tmp(i) = pole(i,k,m)
216 : enddo
217 : #if defined(SPMD)
218 0 : call parcollective(grid%commxy_x,sumop,im,pole_tmp)
219 : #endif
220 0 : pole_sum(k,m) = D0_0
221 0 : do i = 1, im
222 0 : pole_sum(k,m) = pole_sum(k,m) + pole_tmp(i)
223 : enddo
224 : endif
225 : enddo
226 : enddo
227 0 : deallocate( pole_tmp )
228 0 : call t_stopf("pdadj_sumfix")
229 : endif
230 : endif
231 :
232 : ! save results
233 : !$omp parallel do private(i,k,m)
234 7524 : do k = 1, km
235 : ! normalize first
236 2881920 : do m = 1, ntotq+2
237 2881920 : pole_sum(k,m) = pole_sum(k,m)/im
238 : enddo
239 182400 : do i = ifirstxy,ilastxy
240 175104 : delpxy(i,jm,k) = pole_sum(k,1)
241 182400 : ptxy(i,jm,k) = pole_sum(k,2)
242 : enddo
243 2867556 : do m = 1, ntotq
244 71508096 : do i = ifirstxy,ilastxy
245 71500800 : tracer(i,jm,k,m) = pole_sum(k,m+2)
246 : enddo
247 : enddo
248 : enddo
249 :
250 : endif ! jlastxy==jm
251 :
252 : !
253 : ! Fix pe,ps if nontrivial z decomposition
254 : ! Transpose pe - change to better method (16-byte?) later on
255 : !
256 :
257 : !
258 : ! Compute pexy
259 : !
260 : !$omp parallel do private(i, j)
261 58368 : do j = jfirstxy,jlastxy
262 1108992 : do i = ifirstxy, ilastxy
263 1094400 : pexy(i,1,j) = ptop
264 : enddo
265 : enddo
266 :
267 : !$omp parallel do private(i, j, k)
268 58368 : do j = jfirstxy,jlastxy
269 1459200 : do k = 1, km
270 35064576 : do i = ifirstxy, ilastxy
271 35020800 : pexy(i,k+1,j) = pexy(i,k,j) + delpxy(i,j,k)
272 : enddo
273 : enddo
274 : enddo
275 :
276 58368 : do j=jfirstxy,jlastxy
277 1108992 : do i=ifirstxy,ilastxy
278 1094400 : psxy(i,j) = pexy(i,km+1,j)
279 : enddo
280 : enddo
281 :
282 14592 : if (high_alt) then
283 : !$omp parallel do private(i,j,k)
284 0 : do k=2,km
285 0 : do j=jfirstxy,jlastxy
286 0 : do i=ifirstxy,ilastxy
287 0 : cap3vi(i,j,k) = 0.5_r8*(cap3v(i,j,k-1)+cap3v(i,j,k))
288 : enddo
289 : enddo
290 : enddo
291 0 : cap3vi(:,:,1) = 1.5_r8 * cap3v(:,:,1) - 0.5_r8 * cap3v(:,:,2)
292 0 : cap3vi(:,:,km+1) = 1.5_r8 * cap3v(:,:,km) - 0.5_r8 * cap3v(:,:,km-1)
293 : else
294 36611328 : cap3vi(:,:,:) = cap3v(grid%ifirstxy,grid%jfirstxy,1)
295 : endif
296 :
297 : !
298 : ! Update pelnxy and pkxy
299 : !
300 : !$omp parallel do private(i, j, k)
301 58368 : do j=jfirstxy,jlastxy
302 1502976 : do k = 1, km+1
303 36158976 : do i = ifirstxy, ilastxy
304 34670592 : pelnxy(i,k,j) = log( pexy(i,k,j) )
305 34670592 : pkxy (i,j,k) = pexy(i,k,j) ** cap3vi(i,j,k)
306 36115200 : pkln(i,k,j) = log(pkxy(i,j,k))
307 : enddo
308 : enddo
309 : enddo ! jfirstxy:jlastxy loop
310 : !
311 : ! Update pkzxy
312 : !
313 14592 : if (high_alt) then
314 : !$omp parallel do private(i, j, k)
315 0 : do j=jfirstxy,jlastxy
316 0 : do k = 1,km
317 0 : do i = ifirstxy,ilastxy
318 0 : pkzxy(i,j,k) = (pkxy(i,j,k+1)-pkxy(i,j,k))/(pkln(i,k+1,j)-pkln(i,k,j))
319 : enddo
320 : enddo
321 : enddo
322 : else
323 : !$omp parallel do private(i, j, k)
324 58368 : do j=jfirstxy,jlastxy
325 1459200 : do k = 1,km
326 35064576 : do i = ifirstxy,ilastxy
327 35020800 : pkzxy(i,j,k) = (pkxy(i,j,k+1)-pkxy(i,j,k))/(cap3v(i,j,k)*(pelnxy(i,k+1,j)-pelnxy(i,k,j)))
328 : enddo
329 : enddo
330 : enddo
331 : endif
332 :
333 : !
334 : ! Calculate virtual potential temperature
335 : !
336 :
337 : !$omp parallel do private(i, j, k)
338 58368 : do j=jfirstxy,jlastxy
339 1459200 : do k = 1,km
340 35064576 : do i = ifirstxy,ilastxy
341 100859904 : ptxy(i,j,k) = ptxy(i,j,k)* &
342 33619968 : (D1_0+zvir*tracer(i,j,k,1)) &
343 135880704 : /pkzxy(i,j,k)
344 : enddo
345 : enddo
346 : enddo ! jfirstxy:jlastxy loop
347 :
348 : !EOC
349 14592 : end subroutine p_d_adjust
350 : !-----------------------------------------------------------------------
|
