7 use mom_coms,
only : efp_to_real, real_to_efp, efp_list_sum_across_pes
9 use mom_coms,
only : query_efp_overflow_error, reset_efp_overflow_error
15 implicit none ;
private
17 #include <MOM_memory.h>
19 public :: global_i_mean, global_j_mean
20 public :: global_area_mean, global_layer_mean
21 public :: global_area_integral
22 public :: global_volume_mean, global_mass_integral
23 public :: adjust_area_mean_to_zero
28 function global_area_mean(var,G)
30 real,
dimension(SZI_(G), SZJ_(G)),
intent(in) :: var
31 real,
dimension(SZI_(G), SZJ_(G)) :: tmpforsumming
32 real :: global_area_mean
34 integer :: i, j, is, ie, js, je
35 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
37 tmpforsumming(:,:) = 0.
38 do j=js,je ;
do i=is, ie
39 tmpforsumming(i,j) = ( var(i,j) * (g%areaT(i,j) * g%mask2dT(i,j)) )
43 end function global_area_mean
46 function global_area_integral(var,G)
48 real,
dimension(SZI_(G), SZJ_(G)),
intent(in) :: var
49 real,
dimension(SZI_(G), SZJ_(G)) :: tmpforsumming
50 real :: global_area_integral
52 integer :: i, j, is, ie, js, je
53 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
55 tmpforsumming(:,:) = 0.
56 do j=js,je ;
do i=is, ie
57 tmpforsumming(i,j) = ( var(i,j) * (g%areaT(i,j) * g%mask2dT(i,j)) )
61 end function global_area_integral
64 function global_layer_mean(var, h, G, GV)
67 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)),
intent(in) :: var
68 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)),
intent(in) :: h
69 real,
dimension(SZK_(GV)) :: global_layer_mean
71 real,
dimension(SZI_(G), SZJ_(G), SZK_(GV)) :: tmpforsumming, weight
72 real,
dimension(SZK_(GV)) :: scalarij, weightij
73 real,
dimension(SZK_(GV)) :: global_temp_scalar, global_weight_scalar
74 integer :: i, j, k, is, ie, js, je, nz
75 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
77 tmpforsumming(:,:,:) = 0. ; weight(:,:,:) = 0.
79 do k=1,nz ;
do j=js,je ;
do i=is,ie
80 weight(i,j,k) = (gv%H_to_m * h(i,j,k)) * (g%areaT(i,j) * g%mask2dT(i,j))
81 tmpforsumming(i,j,k) = var(i,j,k) * weight(i,j,k)
88 global_layer_mean(k) = scalarij(k) / weightij(k)
91 end function global_layer_mean
94 function global_volume_mean(var, h, G, GV)
97 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
99 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
101 real :: global_volume_mean
104 real,
dimension(SZI_(G), SZJ_(G)) :: tmpforsumming, sum_weight
105 integer :: i, j, k, is, ie, js, je, nz
106 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
108 tmpforsumming(:,:) = 0. ; sum_weight(:,:) = 0.
110 do k=1,nz ;
do j=js,je ;
do i=is,ie
111 weight_here = (gv%H_to_m * h(i,j,k)) * (g%areaT(i,j) * g%mask2dT(i,j))
112 tmpforsumming(i,j) = tmpforsumming(i,j) + var(i,j,k) * weight_here
113 sum_weight(i,j) = sum_weight(i,j) + weight_here
114 enddo ;
enddo ;
enddo
118 end function global_volume_mean
122 function global_mass_integral(h, G, GV, var, on_PE_only)
125 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
127 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
128 optional,
intent(in) :: var
129 logical,
optional,
intent(in) :: on_pe_only
131 real :: global_mass_integral
134 real,
dimension(SZI_(G), SZJ_(G)) :: tmpforsumming
135 logical :: global_sum
136 integer :: i, j, k, is, ie, js, je, nz
137 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
139 tmpforsumming(:,:) = 0.0
141 if (
present(var))
then
142 do k=1,nz ;
do j=js,je ;
do i=is,ie
143 tmpforsumming(i,j) = tmpforsumming(i,j) + var(i,j,k) * &
144 ((gv%H_to_kg_m2 * h(i,j,k)) * (g%areaT(i,j) * g%mask2dT(i,j)))
145 enddo ;
enddo ;
enddo
147 do k=1,nz ;
do j=js,je ;
do i=is,ie
148 tmpforsumming(i,j) = tmpforsumming(i,j) + &
149 ((gv%H_to_kg_m2 * h(i,j,k)) * (g%areaT(i,j) * g%mask2dT(i,j)))
150 enddo ;
enddo ;
enddo
152 global_sum = .true. ;
if (
present(on_pe_only)) global_sum = .not.on_pe_only
156 global_mass_integral = 0.0
157 do j=js,je ;
do i=is,ie
158 global_mass_integral = global_mass_integral + tmpforsumming(i,j)
162 end function global_mass_integral
167 subroutine global_i_mean(array, i_mean, G, mask)
169 real,
dimension(SZI_(G),SZJ_(G)),
intent(in) :: array
170 real,
dimension(SZJ_(G)),
intent(out) :: i_mean
171 real,
dimension(SZI_(G),SZJ_(G)), &
172 optional,
intent(in) :: mask
175 type(
efp_type),
allocatable,
dimension(:) :: asum, mask_sum
177 integer :: is, ie, js, je, idg_off, jdg_off
180 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
181 idg_off = g%idg_offset ; jdg_off = g%jdg_offset
183 call reset_efp_overflow_error()
185 allocate(asum(g%jsg:g%jeg))
186 if (
present(mask))
then
187 allocate(mask_sum(g%jsg:g%jeg))
190 asum(j) = real_to_efp(0.0) ; mask_sum(j) = real_to_efp(0.0)
193 do i=is,ie ;
do j=js,je
194 asum(j+jdg_off) = asum(j+jdg_off) + real_to_efp(array(i,j)*mask(i,j))
195 mask_sum(j+jdg_off) = mask_sum(j+jdg_off) + real_to_efp(mask(i,j))
198 if (query_efp_overflow_error())
call mom_error(fatal, &
199 "global_i_mean overflow error occurred before sums across PEs.")
201 call efp_list_sum_across_pes(asum(g%jsg:g%jeg), g%jeg-g%jsg+1)
202 call efp_list_sum_across_pes(mask_sum(g%jsg:g%jeg), g%jeg-g%jsg+1)
204 if (query_efp_overflow_error())
call mom_error(fatal, &
205 "global_i_mean overflow error occurred during sums across PEs.")
208 mask_sum_r = efp_to_real(mask_sum(j+jdg_off))
209 if (mask_sum_r == 0.0 )
then ; i_mean(j) = 0.0 ;
else
210 i_mean(j) = efp_to_real(asum(j+jdg_off)) / mask_sum_r
216 do j=g%jsg,g%jeg ; asum(j) = real_to_efp(0.0) ;
enddo
218 do i=is,ie ;
do j=js,je
219 asum(j+jdg_off) = asum(j+jdg_off) + real_to_efp(array(i,j))
222 if (query_efp_overflow_error())
call mom_error(fatal, &
223 "global_i_mean overflow error occurred before sum across PEs.")
225 call efp_list_sum_across_pes(asum(g%jsg:g%jeg), g%jeg-g%jsg+1)
227 if (query_efp_overflow_error())
call mom_error(fatal, &
228 "global_i_mean overflow error occurred during sum across PEs.")
231 i_mean(j) = efp_to_real(asum(j+jdg_off)) / real(g%ieg-g%isg+1)
237 end subroutine global_i_mean
241 subroutine global_j_mean(array, j_mean, G, mask)
243 real,
dimension(SZI_(G),SZJ_(G)),
intent(in) :: array
244 real,
dimension(SZI_(G)),
intent(out) :: j_mean
245 real,
dimension(SZI_(G),SZJ_(G)), &
246 optional,
intent(in) :: mask
249 type(
efp_type),
allocatable,
dimension(:) :: asum, mask_sum
251 integer :: is, ie, js, je, idg_off, jdg_off
254 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
255 idg_off = g%idg_offset ; jdg_off = g%jdg_offset
257 call reset_efp_overflow_error()
259 allocate(asum(g%isg:g%ieg))
260 if (
present(mask))
then
261 allocate (mask_sum(g%isg:g%ieg))
264 asum(i) = real_to_efp(0.0) ; mask_sum(i) = real_to_efp(0.0)
267 do i=is,ie ;
do j=js,je
268 asum(i+idg_off) = asum(i+idg_off) + real_to_efp(array(i,j)*mask(i,j))
269 mask_sum(i+idg_off) = mask_sum(i+idg_off) + real_to_efp(mask(i,j))
272 if (query_efp_overflow_error())
call mom_error(fatal, &
273 "global_j_mean overflow error occurred before sums across PEs.")
275 call efp_list_sum_across_pes(asum(g%isg:g%ieg), g%ieg-g%isg+1)
276 call efp_list_sum_across_pes(mask_sum(g%isg:g%ieg), g%ieg-g%isg+1)
278 if (query_efp_overflow_error())
call mom_error(fatal, &
279 "global_j_mean overflow error occurred during sums across PEs.")
282 mask_sum_r = efp_to_real(mask_sum(i+idg_off))
283 if (mask_sum_r == 0.0 )
then ; j_mean(i) = 0.0 ;
else
284 j_mean(i) = efp_to_real(asum(i+idg_off)) / mask_sum_r
290 do i=g%isg,g%ieg ; asum(i) = real_to_efp(0.0) ;
enddo
292 do i=is,ie ;
do j=js,je
293 asum(i+idg_off) = asum(i+idg_off) + real_to_efp(array(i,j))
296 if (query_efp_overflow_error())
call mom_error(fatal, &
297 "global_j_mean overflow error occurred before sum across PEs.")
299 call efp_list_sum_across_pes(asum(g%isg:g%ieg), g%ieg-g%isg+1)
301 if (query_efp_overflow_error())
call mom_error(fatal, &
302 "global_j_mean overflow error occurred during sum across PEs.")
305 j_mean(i) = efp_to_real(asum(i+idg_off)) / real(g%jeg-g%jsg+1)
311 end subroutine global_j_mean
314 subroutine adjust_area_mean_to_zero(array, G, scaling)
316 real,
dimension(SZI_(G),SZJ_(G)),
intent(inout) :: array
317 real,
optional,
intent(out) :: scaling
319 real,
dimension(SZI_(G), SZJ_(G)) :: posvals, negvals, areaxposvals, areaxnegvals
321 real :: areaintposvals, areaintnegvals, posscale, negscale
323 areaxposvals(:,:) = 0.
324 areaxnegvals(:,:) = 0.
326 do j=g%jsc,g%jec ;
do i=g%isc,g%iec
327 posvals(i,j) = max(0., array(i,j))
328 areaxposvals(i,j) = g%areaT(i,j) * posvals(i,j)
329 negvals(i,j) = min(0., array(i,j))
330 areaxnegvals(i,j) = g%areaT(i,j) * negvals(i,j)
336 posscale = 0.0 ; negscale = 0.0
337 if ((areaintposvals>0.).and.(areaintnegvals<0.))
then
338 if (areaintposvals>-areaintnegvals)
then
339 posscale = - areaintnegvals / areaintposvals
340 do j=g%jsc,g%jec ;
do i=g%isc,g%iec
341 array(i,j) = (posscale * posvals(i,j)) + negvals(i,j)
343 elseif (areaintposvals<-areaintnegvals)
then
344 negscale = - areaintposvals / areaintnegvals
345 do j=g%jsc,g%jec ;
do i=g%isc,g%iec
346 array(i,j) = posvals(i,j) + (negscale * negvals(i,j))
350 if (
present(scaling)) scaling = posscale - negscale
352 end subroutine adjust_area_mean_to_zero