FV3DYCORE  Version 2.0.0
fv_fill.F90
Go to the documentation of this file.
1 !***********************************************************************
2 !* GNU Lesser General Public License
3 !*
4 !* This file is part of the FV3 dynamical core.
5 !*
6 !* The FV3 dynamical core is free software: you can redistribute it
7 !* and/or modify it under the terms of the
8 !* GNU Lesser General Public License as published by the
9 !* Free Software Foundation, either version 3 of the License, or
10 !* (at your option) any later version.
11 !*
12 !* The FV3 dynamical core is distributed in the hope that it will be
13 !* useful, but WITHOUT ANYWARRANTY; without even the implied warranty
14 !* of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
15 !* See the GNU General Public License for more details.
16 !*
17 !* You should have received a copy of the GNU Lesser General Public
18 !* License along with the FV3 dynamical core.
19 !* If not, see <http://www.gnu.org/licenses/>.
20 !***********************************************************************
22 !
23 ! Modules Included:
24 ! <table>
25 ! <tr>
26 ! <th>Module Name</th>
27 ! <th>Functions Included</th>
28 ! </tr>
29 ! <tr>
30 ! <td>mpp_domains_mod</td>
31 ! <td>mpp_update_domains, domain2D</td>
32 ! </tr>
33 ! <tr>
34 ! <td>platform_mod</td>
35 ! <td>kind_phys => r8_kind</td>
36 ! </tr>
37 ! </table>
38 
39  use mpp_domains_mod, only: mpp_update_domains, domain2d
40  use platform_mod, only: kind_phys => r8_kind
41 
42  implicit none
43  public fillz
44  public fill_gfs
45  public fill2d
46 
47 contains
48 
51  subroutine fillz(im, km, nq, q, dp)
52  integer, intent(in):: im
53  integer, intent(in):: km
54  integer, intent(in):: nq
55  real , intent(in):: dp(im,km)
56  real , intent(inout) :: q(im,km,nq)
57 ! LOCAL VARIABLES:
58  logical:: zfix(im)
59  real :: dm(km)
60  integer i, k, ic, k1
61  real qup, qly, dup, dq, sum0, sum1, fac
62 
63  do ic=1,nq
64 #ifdef DEV_GFS_PHYS
65 ! Bottom up:
66  do k=km,2,-1
67  k1 = k-1
68  do i=1,im
69  if( q(i,k,ic) < 0. ) then
70  q(i,k1,ic) = q(i,k1,ic) + q(i,k,ic)*dp(i,k)/dp(i,k1)
71  q(i,k ,ic) = 0.
72  endif
73  enddo
74  enddo
75 ! Top down:
76  do k=1,km-1
77  k1 = k+1
78  do i=1,im
79  if( q(i,k,ic) < 0. ) then
80  q(i,k1,ic) = q(i,k1,ic) + q(i,k,ic)*dp(i,k)/dp(i,k1)
81  q(i,k ,ic) = 0.
82  endif
83  enddo
84  enddo
85 #else
86 ! Top layer
87  do i=1,im
88  if( q(i,1,ic) < 0. ) then
89  q(i,2,ic) = q(i,2,ic) + q(i,1,ic)*dp(i,1)/dp(i,2)
90  q(i,1,ic) = 0.
91  endif
92  enddo
93 
94 ! Interior
95  zfix(:) = .false.
96  do k=2,km-1
97  do i=1,im
98  if( q(i,k,ic) < 0. ) then
99  zfix(i) = .true.
100  if ( q(i,k-1,ic) > 0. ) then
101 ! Borrow from above
102  dq = min( q(i,k-1,ic)*dp(i,k-1), -q(i,k,ic)*dp(i,k) )
103  q(i,k-1,ic) = q(i,k-1,ic) - dq/dp(i,k-1)
104  q(i,k ,ic) = q(i,k ,ic) + dq/dp(i,k )
105  endif
106  if ( q(i,k,ic)<0.0 .and. q(i,k+1,ic)>0. ) then
107 ! Borrow from below:
108  dq = min( q(i,k+1,ic)*dp(i,k+1), -q(i,k,ic)*dp(i,k) )
109  q(i,k+1,ic) = q(i,k+1,ic) - dq/dp(i,k+1)
110  q(i,k ,ic) = q(i,k ,ic) + dq/dp(i,k )
111  endif
112  endif
113  enddo
114  enddo
115 
116 ! Bottom layer
117  k = km
118  do i=1,im
119  if( q(i,k,ic)<0. .and. q(i,k-1,ic)>0.) then
120  zfix(i) = .true.
121 ! Borrow from above
122  qup = q(i,k-1,ic)*dp(i,k-1)
123  qly = -q(i,k ,ic)*dp(i,k )
124  dup = min(qly, qup)
125  q(i,k-1,ic) = q(i,k-1,ic) - dup/dp(i,k-1)
126  q(i,k, ic) = q(i,k, ic) + dup/dp(i,k )
127  endif
128  enddo
129 
130 ! Perform final check and non-local fix if needed
131  do i=1,im
132  if ( zfix(i) ) then
133 
134  sum0 = 0.
135  do k=2,km
136  dm(k) = q(i,k,ic)*dp(i,k)
137  sum0 = sum0 + dm(k)
138  enddo
139 
140  if ( sum0 > 0. ) then
141  sum1 = 0.
142  do k=2,km
143  sum1 = sum1 + max(0., dm(k))
144  enddo
145  fac = sum0 / sum1
146  do k=2,km
147  q(i,k,ic) = max(0., fac*dm(k)/dp(i,k))
148  enddo
149  endif
150 
151  endif
152  enddo
153 #endif
154 
155  enddo
156  end subroutine fillz
157 
162  subroutine fill_gfs(im, km, pe2, q, q_min)
163 !SJL: this routine is the equivalent of fillz except that the vertical index is upside down
164  integer, intent(in):: im, km
165  real(kind=kind_phys), intent(in):: pe2(im,km+1)
166  real(kind=kind_phys), intent(in):: q_min
167  real(kind=kind_phys), intent(inout):: q(im,km)
168 ! LOCAL VARIABLES:
169  real(kind=kind_phys) :: dp(im,km)
170  integer:: i, k, k1
171 
172  do k=1,km
173  do i=1,im
174  dp(i,k) = pe2(i,k) - pe2(i,k+1)
175  enddo
176  enddo
177 
178 ! From bottom up:
179  do k=1,km-1
180  k1 = k+1
181  do i=1,im
182  if ( q(i,k)<q_min ) then
183 ! Take mass from above so that q >= q_min
184  q(i,k1) = q(i,k1) + (q(i,k)-q_min)*dp(i,k)/dp(i,k1)
185  q(i,k ) = q_min
186  endif
187  enddo
188  enddo
189 
190 ! From top down:
191  do k=km,2,-1
192  k1 = k-1
193  do i=1,im
194  if ( q(i,k)<0.0 ) then
195 ! Take mass from below
196  q(i,k1) = q(i,k1) + q(i,k)*dp(i,k)/dp(i,k1)
197  q(i,k ) = 0.
198  endif
199  enddo
200  enddo
201 
202  end subroutine fill_gfs
203 
206  subroutine fill2d(is, ie, js, je, ng, km, q, delp, area, domain, bounded_domain, npx, npy)
207 ! This is a diffusive type filling algorithm
208  type(domain2d), intent(INOUT) :: domain
209  integer, intent(in):: is, ie, js, je, ng, km, npx, npy
210  logical, intent(IN):: bounded_domain
211  real, intent(in):: area(is-ng:ie+ng, js-ng:je+ng)
212  real, intent(in):: delp(is-ng:ie+ng, js-ng:je+ng, km)
213  real, intent(inout):: q(is-ng:ie+ng, js-ng:je+ng, km)
214 ! LOCAL VARIABLES:
215  real, dimension(is-ng:ie+ng, js-ng:je+ng,km):: qt
216  real, dimension(is:ie+1, js:je):: fx
217  real, dimension(is:ie, js:je+1):: fy
218  real, parameter:: dif = 0.25
219  integer:: i, j, k
220  integer :: is1, ie1, js1, je1
221 
222  if (bounded_domain) then
223  if (is == 1) then
224  is1 = is-1
225  else
226  is1 = is
227  endif
228  if (ie == npx-1) then
229  ie1 = ie+1
230  else
231  ie1 = ie
232  endif
233  if (js == 1) then
234  js1 = js-1
235  else
236  js1 = js
237  endif
238  if (je == npy-1) then
239  je1 = je+1
240  else
241  je1 = je
242  endif
243  else
244  is1 = is
245  ie1 = ie
246  js1 = js
247  je1 = je
248  endif
249 
250 !$OMP parallel do default(shared)
251  do k=1, km
252  do j=js1,je1
253  do i=is1,ie1
254  qt(i,j,k) = q(i,j,k)*delp(i,j,k)*area(i,j)
255  enddo
256  enddo
257  enddo
258  call mpp_update_domains(qt, domain, whalo=1, ehalo=1, shalo=1, nhalo=1)
259 
260 !$OMP parallel do default(shared) private(fx,fy)
261  do k=1, km
262  fx(:,:) = 0.
263  do j=js,je
264  do i=is,ie+1
265  if( qt(i-1,j,k)*qt(i,j,k)<0. ) fx(i,j) = qt(i-1,j,k) - qt(i,j,k)
266  enddo
267  enddo
268  fy(:,:) = 0.
269  do j=js,je+1
270  do i=is,ie
271  if( qt(i,j-1,k)*qt(i,j,k)<0. ) fy(i,j) = qt(i,j-1,k) - qt(i,j,k)
272  enddo
273  enddo
274  do j=js,je
275  do i=is,ie
276  q(i,j,k) = q(i,j,k)+dif*(fx(i,j)-fx(i+1,j)+fy(i,j)-fy(i,j+1))/(delp(i,j,k)*area(i,j))
277  enddo
278  enddo
279  enddo
280 
281  end subroutine fill2d
282 
283 end module fv_fill_mod
subroutine, public fillz(im, km, nq, q, dp)
The subroutine &#39;fillz&#39; is for mass-conservative filling of nonphysical negative values in the tracers...
Definition: fv_fill.F90:52
subroutine, public fill_gfs(im, km, pe2, q, q_min)
The subroutine &#39;fill_gfs&#39; is for mass-conservative filling of nonphysical negative values in the trac...
Definition: fv_fill.F90:163
subroutine, public fill2d(is, ie, js, je, ng, km, q, delp, area, domain, bounded_domain, npx, npy)
The subroutine &#39;fill2D&#39; fills in nonphysical negative values in a single scalar field using a two-dim...
Definition: fv_fill.F90:207