FMS  2025.04
Flexible Modeling System
fms_diag_time_reduction.F90
1 !***********************************************************************
2 !* Apache License 2.0
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* Licensed under the Apache License, Version 2.0 (the "License");
7 !* you may not use this file except in compliance with the License.
8 !* You may obtain a copy of the License at
9 !*
10 !* http://www.apache.org/licenses/LICENSE-2.0
11 !*
12 !* FMS is distributed in the hope that it will be useful, but WITHOUT
13 !* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied;
14 !* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
15 !* PARTICULAR PURPOSE. See the License for the specific language
16 !* governing permissions and limitations under the License.
17 !***********************************************************************
18 
19 !> @defgroup fms_diag_time_reduction_mod fms_diag_time_reduction_mod
20 !> @ingroup diag_manager
21 !> @brief fms_diag_time_reduction_mod defines classes encapsulating the diag_manager
22 !! time redution types.
23 !!
24 !> @author Miguel Zuniga
25 !!
26 !> @file
27 !> @brief File for @ref fms_diag_time_reduction_mod
28 !> @addtogroup fms_diag_time_reduction_mod
29 !> @{
30 MODULE fms_diag_time_reduction_mod
31 
32  USE diag_data_mod, only: every_time
33  USE fms_mod, ONLY: error_mesg, fatal
34 
35  implicit none
36 
37  !!TODO: (Future effort) Note that time_diurnal processing is a little different
38  !! and more complex than the other reduction methods, and therefore refactoring its
39  !! processing may simplify the overall related codebase. The refactoring,
40  !! if possible, may be done elsewhere in the diag_manager.
41 
42  !!These parametes are the possible kinds of time reduction operations.
43  INTEGER, PARAMETER :: time_none = 0 !< There is no reduction method
44  INTEGER, PARAMETER :: time_average = 1 !< The reduction method is average
45  INTEGER, PARAMETER :: time_rms = 2 !< The reduction method is rms
46  INTEGER, PARAMETER :: time_max = 3 !< The reduction method is max
47  INTEGER, PARAMETER :: time_min = 4 !< The reduction method is min
48  INTEGER, PARAMETER :: time_sum = 5 !< The reudction method is sum
49  INTEGER, PARAMETER :: time_diurnal = 6 !< The reduction method is diurnal
50  INTEGER, PARAMETER :: time_power = 7 !< The reduction method is power
51 
52  !> @brief Class fmsDiagTimeReduction_type has an encapsulation of the "Fortran enum" time
53  !! reduction integer parameters, plus an encapsulation of the groupings of
54  !! the time reduction types. It is intended to provide some of the functionality
55  !! that was coded in the legacy function diag_data.F90:init_output_fields.
56  !! The functionality in the end is used by send_data in (EFFICIENT) do loops calling
57  !! the weighting or math functions to update buffers.
58  !! The integer parameters above are the legal time reduction types,
59  !! but they are not necessarily mutually exclusive in some contexts.
60  !!
61  !> @addtogroup fms_diag_time_reduction_mod
63  integer , private :: the_time_reduction !< The time reduction type, as an integer defined above.
64  logical , private :: time_averaging !< Set true iff time_average, time_rms, time_power or time_diurnal is true
65  logical , private :: time_ops !< Set true iff time_min, time_max, time_rms or time_average is true.
66  CONTAINS
67  procedure, public :: do_time_averaging => do_time_averaging_imp
68  procedure, public :: has_time_ops => has_time_ops_imp
69  procedure, public :: is_time_none => is_time_none_imp
70  procedure, public :: is_time_average => is_time_average_imp
71  procedure, public :: is_time_rms => is_time_rms_imp
72  procedure, public :: is_time_max => is_time_max_imp
73  procedure, public :: is_time_min => is_time_min_imp
74  procedure, public :: is_time_sum => is_time_sum_imp
75  procedure, public :: is_time_diurnal => is_time_diurnal_imp
76  procedure, public :: is_time_power => is_time_power_imp
77  procedure, public :: initialize
78  procedure, public :: copy
80 
81  !> @brief This interface is for the class constructor.
82  !> @addtogroup fms_diag_time_reduction_mod
85  end interface fmsdiagtimereduction_type
86 
87 CONTAINS
88 
89  !> @brief The class contructors. Just allocates the class and calls an initializer
90  !! @return An allocated instance of fmsDiagTimeReduction_type, which is nitialized using
91  !! provided values for arguments dt and out_freqeuncy.
92  function fmsdiagtimereduction_type_constructor(dt, out_frequency) result(time_redux)
93  integer, intent(in) :: dt !< The redution type (time_rms, time_power, etc)
94  integer, intent(in) :: out_frequency !< The output frequency.
95  class(fmsdiagtimereduction_type), allocatable :: time_redux !< The instance of the fmsDiagTimeReduction_type
96  !!class allocated and returned by this constructor.
97  allocate(time_redux)
98  call time_redux%initialize(dt, out_frequency)
100 
101  !> @brief Initialize the object. As an alternative to the constructor, one can
102  !! allocate an fmsDiagTimeReduction_type instance, then call its initialize function.
103  subroutine initialize(this, dt, out_frequency)
104  class(fmsdiagtimereduction_type), intent(inout) :: this !< The fmsDiagTimeReduction_type object
105  integer, intent(in) :: dt !< The redution type (time_rms, time_porer, etc)
106  integer, intent(in) :: out_frequency !< The output frequency.
107 
108  this%the_time_reduction = dt
109 
110  !! Set the time_averaging flag
111  !! See legacy init_ouput_fields function, lines 1470ff
112  IF(( dt .EQ. time_average) .OR. (dt .EQ. time_rms) .OR. (dt .EQ. time_power) .OR. &
113  & (dt .EQ. time_diurnal)) THEN
114  this%time_averaging = .true.
115  ELSE
116  this%time_averaging= .false.
117  IF((dt .NE. time_max) .AND. (dt .ne. time_min) .AND. (dt .NE. time_sum) &
118  & .AND. (dt .NE. time_none)) THEN
119  CALL error_mesg('fmsDiagTimeReduction_type: initialize', &
120  & 'time_averaging=.false. but reduction type not compatible', fatal)
121  ENDIF
122  END IF
123 
124  !!TODO: (MDM) Add other checks? E.g. If time_averaging == .false., then
125  !! out_frequency == EVERY_TIME
126 
127  IF((dt .EQ. time_min) .OR. (dt .EQ. time_max) .OR. &
128  & ( dt .EQ. time_average) .OR. (dt .EQ. time_sum) ) THEN
129  this%time_ops = .true.
130  ELSE
131  this%time_ops = .false.
132  END IF
133  end subroutine initialize
134 
135  !> @brief Copy the source time reduction object into the this object.
136  subroutine copy(this, source)
137  class(fmsdiagtimereduction_type),intent(inout):: this !< The fmsDiagTimeReduction_type object
138  class(fmsdiagtimereduction_type),intent(in):: source !< The fmsDiagTimeReduction_type object
139  this%the_time_reduction = source%the_time_reduction
140  this%time_averaging = source%time_averaging
141  this%time_ops = source%time_ops
142  end subroutine copy
143 
144  !> \brief Returns true if any of time_min, time_max, time_rms or time_average is true.
145  !! @return true if any of time_min, time_max, time_rms or time_average is true.
146  pure function has_time_ops_imp (this)
147  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
148  logical :: has_time_ops_imp
149  has_time_ops_imp = this%time_ops
150  end function has_time_ops_imp
151 
152  !> \brief Returns true iff time_averaging is true.
153  !! @return true iff time_averaging is true.
154  pure function do_time_averaging_imp (this)
155  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
156  logical :: do_time_averaging_imp
157  do_time_averaging_imp = this%time_averaging
158  end function do_time_averaging_imp
159 
160  !> \brief Returns true iff the_time_reduction is time_average
161  !! @return true iff the_time_reduction is time_average
162  pure function is_time_average_imp (this)
163  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
164  logical :: is_time_average_imp
165  is_time_average_imp = this%the_time_reduction .EQ. time_average
166  end function is_time_average_imp
167 
168  !> \brief Returns true iff the_time_reduction is time_none
169  !! @return true iff the_time_reduction is time_none
170  pure function is_time_none_imp (this)
171  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
172  logical :: is_time_none_imp
173  is_time_none_imp = (this%the_time_reduction .EQ. time_none)
174  end function is_time_none_imp
175 
176  !> \brief Returns true iff the_time_reduction is time_rms
177  !! @return true iff the_time_reduction is time_rms
178  pure function is_time_rms_imp (this)
179  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
180  logical :: is_time_rms_imp
181  is_time_rms_imp = this%the_time_reduction .EQ. time_rms
182  end function is_time_rms_imp
183 
184  !> \brief Returns true iff the_time_reduction is time_max
185  !! @return true iff the_time_reduction is time_max
186  pure function is_time_max_imp (this)
187  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
188  logical :: is_time_max_imp
189  is_time_max_imp = this%the_time_reduction .EQ. time_max
190  end function is_time_max_imp
191 
192  !> \brief Returns true iff the_time_reduction is time_min
193  !! @return true iff the_time_reduction is time_min
194  pure function is_time_min_imp (this)
195  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
196  logical :: is_time_min_imp
197  is_time_min_imp = this%the_time_reduction .EQ. time_min
198  end function is_time_min_imp
199 
200  !> \brief Returns true iff the_time_reduction is time_sum
201  !! @return true iff the_time_reduction is time_sum
202  pure function is_time_sum_imp (this)
203  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
204  logical :: is_time_sum_imp
205  is_time_sum_imp = this%the_time_reduction .EQ. time_sum
206  end function is_time_sum_imp
207 
208  !> \brief Returns true iff the_time_reduction is time_diurnal
209  !! @return true iff the_time_reduction is time_diurnal
210  pure function is_time_diurnal_imp (this)
211  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
212  logical :: is_time_diurnal_imp
213  is_time_diurnal_imp = this%the_time_reduction .EQ. time_diurnal
214  end function is_time_diurnal_imp
215 
216  !> \brief Returns true iff the_time_reduction is time_power
217  !! @return true iff the_time_reduction is time_power
218  pure function is_time_power_imp (this)
219  class(fmsdiagtimereduction_type), intent(in) :: this !<The object this function is bound to.
220  logical :: is_time_power_imp
221  is_time_power_imp = this%the_time_reduction .EQ. time_power
222  end function is_time_power_imp
223 
224 end module fms_diag_time_reduction_mod
225 !> @}
226 ! close documentation grouping
pure logical function is_time_none_imp(this)
Returns true iff the_time_reduction is time_none.
pure logical function is_time_sum_imp(this)
Returns true iff the_time_reduction is time_sum.
pure logical function do_time_averaging_imp(this)
Returns true iff time_averaging is true.
pure logical function is_time_min_imp(this)
Returns true iff the_time_reduction is time_min.
class(fmsdiagtimereduction_type) function, allocatable fmsdiagtimereduction_type_constructor(dt, out_frequency)
The class contructors. Just allocates the class and calls an initializer.
pure logical function is_time_rms_imp(this)
Returns true iff the_time_reduction is time_rms.
pure logical function is_time_max_imp(this)
Returns true iff the_time_reduction is time_max.
subroutine copy(this, source)
Copy the source time reduction object into the this object.
integer, parameter time_min
The reduction method is min.
integer, parameter time_diurnal
The reduction method is diurnal.
integer, parameter time_power
The reduction method is power.
pure logical function is_time_power_imp(this)
Returns true iff the_time_reduction is time_power.
pure logical function is_time_diurnal_imp(this)
Returns true iff the_time_reduction is time_diurnal.
integer, parameter time_average
The reduction method is average.
integer, parameter time_sum
The reudction method is sum.
integer, parameter time_rms
The reduction method is rms.
subroutine initialize(this, dt, out_frequency)
Initialize the object. As an alternative to the constructor, one can allocate an fmsDiagTimeReduction...
pure logical function is_time_average_imp(this)
Returns true iff the_time_reduction is time_average.
integer, parameter time_max
The reduction method is max.
pure logical function has_time_ops_imp(this)
Returns true if any of time_min, time_max, time_rms or time_average is true.
subroutine, public error_mesg(routine, message, level)
Print notes, warnings and error messages; terminates program for warning and error messages....
Definition: fms.F90:441