FMS 2025.01-dev
Flexible Modeling System
Loading...
Searching...
No Matches
diag_integral.inc
1!***********************************************************************
2!* GNU Lesser General Public License
3!*
4!* This file is part of the GFDL Flexible Modeling System (FMS).
5!*
6!* FMS is free software: you can redistribute it and/or modify it under
7!* the terms of the GNU Lesser General Public License as published by
8!* the Free Software Foundation, either version 3 of the License, or (at
9!* your option) any later version.
10!*
11!* FMS is distributed in the hope that it will be useful, but WITHOUT
12!* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13!* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14!* for more details.
15!*
16!* You should have received a copy of the GNU Lesser General Public
17!* License along with FMS. If not, see <http://www.gnu.org/licenses/>.
18!***********************************************************************
19!> @defgroup diag_integral_mod diag_integral_mod
20!> @ingroup diag_integral
21!!
22!! @author Fei Liu <Fei.Liu@noaa.gov>
23!!
24!! @brief This module computes and outputs global and / or hemispheric physics
25!! integrals.
26!###############################################################################
27!> @brief Perform a 2 dimensional summation of named field
28!!
29!! @implements sum_diag_integral_field
30!!
31!! <b> Template: </b>
32!!
33!! @code{.f90}
34!! call sum_field_2d (name, field_data2d, is, js)
35!! @endcode
36!!
37!! <b> Parameters: </b>
38!!
39!! @code{.f90}
40!! character(len=*), intent(in) :: name
41!! real, intent(in) :: field_data2d(:,:)
42!! integer, optional, intent(in) :: is, js
43!! @endcode
44!!
45!! @param [in] <name> Name of the field to be integrated
46!! @param [in] <field_data2d> field of integrands to be summed over
47!! @param [in] <is, js> starting i,j indices over which summation is to occur
48!!
49subroutine sum_field_2d_ (name, field_data2d, is, js)
50
51character(len=*), intent(in) :: name !< Name of the field to be integrated
52real(FMS_DI_KIND_), intent(in) :: field_data2d(:,:) !< field of integrands to be summed over
53integer, optional, intent(in) :: is !< starting i indices over which summation is to occur
54integer, optional, intent(in) :: js !< starting j indices over which summation is to occur
55
56!-------------------------------------------------------------------------------
57! local variables:
58!-------------------------------------------------------------------------------
59 integer :: field !< index of desired integral
60 integer :: i1 !< location indices of current data in
61 !! processor-global coordinates
62 integer :: j1 !< location indices of current data in
63 !! processor-global coordinates
64 integer :: i2 !< location indices of current data in
65 !! processor-global coordinates
66 integer :: j2 !< location indices of current data in
67 !! processor-global coordinates
68
69
70!-------------------------------------------------------------------------------
71! be sure module has been initialized.
72!-------------------------------------------------------------------------------
73 if (.not. module_is_initialized ) then
74 call error_mesg ('diag_integral_mod', &
75 'module has not been initialized', fatal )
76 endif
77
78!-------------------------------------------------------------------------------
79! obtain the index of the current integral. make certain it is valid.
80!-------------------------------------------------------------------------------
81 field = get_field_index(name)
82 if (field == 0) then
83 call error_mesg ('diag_integral_mod', &
84 'field does not exist', fatal)
85 endif
86
87!-------------------------------------------------------------------------------
88! define the processor-global indices of the current data. use the
89! value 1 for the initial grid points, if is and js are not input.
90!-------------------------------------------------------------------------------
91 i1 = 1; if (present(is)) i1 = is
92 j1 = 1; if (present(js)) j1 = js
93 i2 = i1 + size(field_data2d,1) - 1
94 j2 = j1 + size(field_data2d,2) - 1
95
96!-------------------------------------------------------------------------------
97! increment the count of points toward this integral and add the
98! values at this set of grid points to the accumulation array.
99!-------------------------------------------------------------------------------
100!$OMP CRITICAL
101 field_count(field) = field_count(field) + &
102 size(field_data2d,1)*size(field_data2d,2)
103 field_sum(field) = field_sum(field) + &
104 sum(real(field_data2d,r8_kind) * area(i1:i2,j1:j2))
105
106!$OMP END CRITICAL
107
108end subroutine sum_field_2d_
109
110
111
112!###############################################################################
113!> @brief Perform a 3 dimensional summation of named field
114!!
115!! @implements sum_diag_integral_field
116!!
117!! <b> Template: </b>
118!!
119!! @code{.f90}
120!! call sum_field_3d (name, field_data3d, is, js)
121!! @endcode
122!!
123!! <b> Parameters: </b>
124!!
125!! @code{.f90}
126!! character(len=*), intent(in) :: name
127!! real, intent(in) :: field_data3d(:,:,:)
128!! integer, optional, intent(in) :: is, js
129!! @endcode
130!!
131!! @param [in] <name> Name of the field to be integrated
132!! @param [in] <field_data3d> field of integrands to be summed over
133!! @param [in] <is, js> starting i,j indices over which summation is to occur
134!!
135subroutine sum_field_3d_ (name, field_data3d, is, js)
136
137character(len=*), intent(in) :: name !< Name of the field to be integrated
138real(FMS_DI_KIND_), intent(in) :: field_data3d(:,:,:) !< field of integrands to be summed over
139integer, optional, intent(in) :: is !< starting i,j indices over which summation is to occur
140integer, optional, intent(in) :: js !< starting i,j indices over which summation is to occur
141
142!-------------------------------------------------------------------------------
143! local variables:
144! data2
145! field ! index of desired integral
146! i1, j1, i2, j2 ! location indices of current data in
147! processor-global coordinates
148!-------------------------------------------------------------------------------
149 real(r8_kind), dimension (size(field_data3d,1), & size(field_data3d,2)) :: data2
150
151 integer :: field !< index of desired integral
152 integer :: i1 !< location indices of current data in
153 !! processor-global coordinates
154 integer :: j1 !< location indices of current data in
155 !! processor-global coordinates
156 integer :: i2 !< location indices of current data in
157 !! processor-global coordinates
158 integer :: j2 !< location indices of current data in
159 !! processor-global coordinates
160
161
162!-------------------------------------------------------------------------------
163! be sure module has been initialized.
164!-------------------------------------------------------------------------------
165 if (.not. module_is_initialized ) then
166 call error_mesg ('diag_integral_mod', &
167 'module has not been initialized', fatal )
168 endif
169
170!-------------------------------------------------------------------------------
171! obtain the index of the current integral. make certain it is valid.
172!-------------------------------------------------------------------------------
173 field = get_field_index(name)
174 if (field == 0) then
175 call error_mesg ('diag_integral_mod', &
176 'field does not exist', fatal)
177 endif
178
179!-------------------------------------------------------------------------------
180! define the processor-global indices of the current data. use the
181! value 1 for the initial grid points, if is and js are not input.
182!-------------------------------------------------------------------------------
183 i1 = 1; if (present(is)) i1 = is
184 j1 = 1; if (present(js)) j1 = js
185 i2 = i1 + size(field_data3d,1) - 1
186 j2 = j1 + size(field_data3d,2) - 1
187
188!-------------------------------------------------------------------------------
189! increment the count of points toward this integral. sum first
190! in the vertical and then add the values at this set of grid points
191! to the accumulation array.
192!-------------------------------------------------------------------------------
193!$OMP CRITICAL
194 field_count(field) = field_count(field) + &
195 size(field_data3d,1)*size(field_data3d,2)
196 data2 = sum(real(field_data3d,r8_kind),3)
197 field_sum(field) = field_sum(field) + &
198 sum(data2 * area(i1:i2,j1:j2))
199
200!$OMP END CRITICAL
201
202end subroutine sum_field_3d_
203
204
205
206!###############################################################################
207!> @brief Perform a 3 dimensional weighted summation of named field
208!!
209!! @implements sum_diag_integral_field
210!!
211!! <b> Template: </b>
212!!
213!! @code{.f90}
214!! call sum_field_wght_3d (name, field_data3d, wt, is, js)
215!! @endcode
216!!
217!! <b> Parameters: </b>
218!!
219!! @code{.f90}
220!! character(len=*), intent(in) :: name
221!! real, intent(in) :: field_data3d(:,:,:), wt(:,:,:)
222!! integer, optional, intent(in) :: is, js
223!! @endcode
224!!
225!! @param [in] <name> Name of the field to be integrated
226!! @param [in] <field_data3d> field of integrands to be summed over
227!! @param [in] <wt> the weight function to be evaluated at summation
228!! @param [in] <is, js> starting i,j indices over which summation is to occur
229!!
230subroutine sum_field_wght_3d_ (name, field_data3d, wt, is, js)
231
232character(len=*), intent(in) :: name !< Name of the field to be integrated
233real(FMS_DI_KIND_), intent(in) :: field_data3d(:,:,:) !< field of integrands to be summed over
234real(FMS_DI_KIND_), intent(in) :: wt(:,:,:) !< the weight function to be evaluated at summation
235integer, optional, intent(in) :: is !< starting i indices over which summation is to occur
236integer, optional, intent(in) :: js !< starting j indices over which summation is to occur
237
238!-------------------------------------------------------------------------------
239! local variables:
240! data2
241! field ! index of desired integral
242! i1, j1, i2, j2 ! location indices of current data in
243! processor-global coordinates
244!-------------------------------------------------------------------------------
245 real(r8_kind), dimension (size(field_data3d,1),size(field_data3d,2)) :: data2
246 integer :: field !< index of desired integral
247 integer :: i1 !< location indices of current data in
248 !! processor-global coordinates
249 integer :: j1 !< location indices of current data in
250 !! processor-global coordinates
251 integer :: i2 !< location indices of current data in
252 !! processor-global coordinates
253 integer :: j2 !< location indices of current data in
254 !! processor-global coordinates
255
256!-------------------------------------------------------------------------------
257! be sure module has been initialized.
258!-------------------------------------------------------------------------------
259 if (.not. module_is_initialized ) then
260 call error_mesg ('diag_integral_mod', &
261 'module has not been initialized', fatal )
262 endif
263
264!-------------------------------------------------------------------------------
265! obtain the index of the current integral. make certain it is valid.
266!-------------------------------------------------------------------------------
267 field = get_field_index(name)
268 if (field == 0) then
269 call error_mesg ('diag_integral_mod', &
270 'field does not exist', fatal)
271 endif
272
273!-------------------------------------------------------------------------------
274! define the processor-global indices of the current data. use the
275! value 1 for the initial grid points, if is and js are not input.
276!-------------------------------------------------------------------------------
277 i1 = 1; if (present(is)) i1 = is
278 j1 = 1; if (present(js)) j1 = js
279 i2 = i1 + size(field_data3d,1) - 1
280 j2 = j1 + size(field_data3d,2) - 1
281
282!-------------------------------------------------------------------------------
283! increment the count of points toward this integral. sum first
284! in the vertical (including a vertical weighting factor) and then
285! add the values at this set of grid points to the accumulation
286! array.
287!-------------------------------------------------------------------------------
288!$OMP CRITICAL
289 field_count(field) = field_count(field) + &
290 size(field_data3d,1)*size(field_data3d,2)
291 data2 = vert_diag_integral(real(field_data3d,r8_kind), real(wt,r8_kind))
292 field_sum(field) = field_sum(field) + &
293 sum(data2 * area(i1:i2,j1:j2))
294
295!$OMP END CRITICAL
296
297end subroutine sum_field_wght_3d_
298
299
300
301!###############################################################################
302!> @brief Perform a 2 dimensional hemispherical summation of named field
303!!
304!! @implements sum_diag_integral_field
305!!
306!! <b> Template: </b>
307!!
308!! @code{.f90}
309!! call sum_field_2d_hemi (name, field_data2d, is, ie, js, je)
310!! @endcode
311!!
312!! <b> Parameters: </b>
313!!
314!! @code{.f90}
315!! character(len=*), intent(in) :: name
316!! real, intent(in) :: field_data2d(:,:)
317!! integer, intent(in) :: is, js, ie, je
318!! @endcode
319!!
320!! @param [in] <name> Name of the field to be integrated
321!! @param [in] <field_data2d> field of integrands to be summed over
322!! @param [in] <is, js, ie, je> starting/ending i,j indices over which summation
323!! is to occur
324!!
325subroutine sum_field_2d_hemi_ (name, field_data2d, is, ie, js, je)
326
327character(len=*), intent(in) :: name !< Name of the field to be integrated
328real(FMS_DI_KIND_),intent(in) :: field_data2d(:,:) !< field of integrands to be summed over
329integer, intent(in) :: is !< starting/ending i,j indices over which summation
330 !! is to occur
331integer, intent(in) :: js !< starting/ending i,j indices over which summation
332 !! is to occur
333integer, intent(in) :: ie !< starting/ending i,j indices over which summation
334 !! is to occur
335integer, intent(in) :: je !< starting/ending i,j indices over which summation
336 !! is to occur
337
338!-------------------------------------------------------------------------------
339! local variables:
340! field ! index of desired integral
341! i1, j1, i2, j2 ! location indices of current data in
342! processor-global coordinates
343!-------------------------------------------------------------------------------
344 integer :: field !< index of desired integral
345 integer :: i1 !< location indices of current data in
346 !! processor-global coordinates
347 integer :: j1 !< location indices of current data in
348 !! processor-global coordinates
349 integer :: i2 !< location indices of current data in
350 !! processor-global coordinates
351 integer :: j2 !< location indices of current data in
352 !! processor-global coordinates
353
354!-------------------------------------------------------------------------------
355! be sure module has been initialized.
356!-------------------------------------------------------------------------------
357 if (.not. module_is_initialized ) then
358 call error_mesg ('diag_integral_mod', &
359 'module has not been initialized', fatal )
360 endif
361
362!-------------------------------------------------------------------------------
363! obtain the index of the current integral. make certain it is valid.
364!-------------------------------------------------------------------------------
365 field = get_field_index(name)
366 if (field == 0) then
367 call error_mesg ('diag_integral_mod', &
368 'field does not exist', fatal)
369 endif
370
371!-------------------------------------------------------------------------------
372! define the processor-global indices of the current data. this form
373! is needed to handle case of 2d domain decomposition with physics
374! window smaller than processor domain size.
375!-------------------------------------------------------------------------------
376 i1 = mod( (is-1), size(field_data2d,1) ) + 1
377 i2 = i1 + size(field_data2d,1) - 1
378
379!-------------------------------------------------------------------------------
380! for a hemispheric sum, sum one jrow at a time in case a processor
381! has data from both hemispheres.
382!-------------------------------------------------------------------------------
383 j1 = mod( (js-1) ,size(field_data2d,2) ) + 1
384 j2 = j1
385
386!-------------------------------------------------------------------------------
387! increment the count of points toward this integral. include hemi-
388! spheric factor of 2 in field_count. add the data values at this
389! set of grid points to the accumulation array.
390!-------------------------------------------------------------------------------
391!$OMP CRITICAL
392 field_count(field) = field_count(field) + 2* (i2-i1+1)*(j2-j1+1)
393 field_sum(field) = field_sum(field) + &
394 sum( real(field_data2d(i1:i2,j1:j2),r8_kind) * area(is:ie,js:je))
395
396!$OMP END CRITICAL
397
398end subroutine sum_field_2d_hemi_
399
400
401!> @}
402! close documentation grouping
403