1 !***********************************************************************
2 !* GNU Lesser General Public License
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
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.
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
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 field_manager_mod field_manager_mod
20 !> @ingroup field_manager
21 !> @brief Reads entries from a field table and stores this
22 !! information along with the type of field it belongs to.
24 !> This allows the component models to query the field manager to see if non-default
25 !! methods of operation are desired. In essence the field table is a
26 !! powerful type of namelist. Default values can be provided for all the
27 !! fields through a namelist, individual fields can be modified through
28 !! the field table however.
30 !> @author William Cooke
32 !! An example of field table entries could be
34 !! "tracer","atmos_mod","sphum"
36 !! "tracer","atmos_mod","sf6"
37 !! "longname","sulf_hex"
38 !! "advection_scheme_horiz","2nd_order"
39 !! "Profile_type","Fixed","surface_value = 0.0E+00"/
41 !! "prog_tracers","ocean_mod","age_global"
42 !! horizontal-advection-scheme = mdfl_sweby
43 !! vertical-advection-scheme = mdfl_sweby
44 !! restart_file = ocean_age.res.nc
47 !! The field table consists of entries in the following format.
49 !! The first line of an entry should consist of three quoted strings.
51 !! The first quoted string will tell the field manager what type of
54 !! The second quoted string will tell the field manager which model the
55 !! field is being applied to.
56 !! The supported types at present are
58 !! "coupler_mod" for the coupler,
59 !! "atmos_mod" for the atmosphere model,
60 !! "ocean_mod" for the ocean model,
61 !! "land_mod" for the land model, and,
62 !! "ice_mod" for the ice model.
64 !! The third quoted string should be a unique name that can be used as a
67 !! The second and following lines of each entry are called methods in
68 !! this context. Methods can be developed within any module and these
69 !! modules can query the field manager to find any methods that are
70 !! supplied in the field table.
72 !! These lines can be coded quite flexibly.
74 !! The line can consist of two or three quoted strings or a simple unquoted
77 !! If the line consists two or three quoted strings, then the first string will
78 !! be an identifier that the querying module will ask for.
80 !! The second string will be a name that the querying module can use to
81 !! set up values for the module.
83 !! The third string, if present, can supply parameters to the calling module that can be
84 !! parsed and used to further modify values.
86 !! If the line consists of a simple unquoted string then quotes are not allowed
87 !! in any part of the line.
89 !! An entry is ended with a backslash (/) as the final character in a
92 !! Comments can be inserted in the field table by having a # as the
93 !! first character in the line.
95 !! In the example above we have three field entries.
97 !! The first is a simple declaration of a tracer called "sphum".
99 !! The second is for a tracer called "sf6". In this case a field named
100 !! "longname" will be given the value "sulf_hex". A field named
101 !! "advection_scheme_horiz" will be given the value "2nd_order". Finally a field
102 !! name "Profile_type" will be given a child field called "Fixed", and that field
103 !! will be given a field called "surface_value" with a real value of 0.0E+00.
105 !! The third entry is an example of a oceanic age tracer. Note that the
106 !! method lines are formatted differently here. This is the flexibility mentioned
109 !! With these formats, a number of restrictions are required.
111 !! The following formats are equally valid.
113 !! "longname","sulf_hex"
114 !! "longname = sulf_hex"
115 !! longname = sulf_hex
117 !! However the following is not valid.
119 !! longname = "sulf_hex"
122 !! In the SF6 example above the last line of the entry could be written in the
125 !! "Profile_type","Fixed","surface_value = 0.0E+00"/
126 !! Profile_type/Fixed/surface_value = 0.0E+00/
129 !! Values supplied with fields are converted to the various types with the
130 !! following assumptions.
132 !! Real values : These values contain a decimal point or are in exponential format.
133 !! These values only support e or E format for exponentials.
134 !! e.g. 10.0, 1e10 and 1E10 are considered to be real numbers.
136 !! Integer values : These values only contain numbers.
137 !! e.g 10 is an integer. 10.0 and 1e10 are not.
139 !! Logical values : These values are supplied as one of the following formats.
140 !! T, .T., TRUE, .TRUE.
141 !! t, .t., true, .true.
142 !! F, .F., FALSE, .FALSE.
143 !! f, .f., false, .false.
144 !! These will be converted to T or F in a dump of the field.
146 !! Character strings : These values are assumed to be strings if a character
147 !! other than an e (or E) is in the value. Numbers can be suppled in the value.
148 !! If the value does not meet the criteria for a real, integer or logical type,
149 !! it is assumed to be a character type.
151 !! The entries within the field table can be designed by the individual
152 !! authors of code to allow modification of their routines.
155 !> @addtogroup field_manager_mod
157 module field_manager_mod
158 !TODO this variable can be removed when the legacy table is no longer used
160 #define MAXFIELDS_ 250
163 !TODO this variable can be removed when the legacy table is not longer used
164 #ifndef MAXFIELDMETHODS_
165 #define MAXFIELDMETHODS_ 250
169 ! <CONTACT EMAIL="William.Cooke@noaa.gov"> William Cooke
172 ! <REVIEWER EMAIL="Richard.Slater@noaa.gov"> Richard D. Slater
175 ! <REVIEWER EMAIL="Matthew.Harrison@noaa.gov"> Matthew Harrison
178 ! <REVIEWER EMAIL="John.Dunne@noaa.gov"> John P. Dunne
181 use mpp_mod, only : mpp_error, &
190 use fms_mod, only : lowercase, &
191 write_version_number, &
193 use fms2_io_mod, only: file_exists
194 use platform_mod, only: r4_kind, r8_kind
202 #include<file_version.h>
203 logical :: module_is_initialized = .false.
205 public :: field_manager_init !< (nfields, [table_name]) returns number of fields
206 public :: field_manager_end !< ()
207 public :: find_field_index !< (model, field_name) or (list_path)
208 public :: find_field_index_old !< (model, field_name) returns index of field_name in
209 public :: find_field_index_new
210 public :: get_field_info !< (n,fld_type,fld_name,model,num_methods)
211 !! Returns parameters relating to field n.
212 public :: get_field_method !< (n, m, method) Returns the m-th method of field n
213 public :: get_field_methods !< (n, methods) Returns the methods related to field n
214 public :: parse !< (text, label, values) Overloaded function to parse integer,
215 !! real or character. Parse returns the number of values
216 !! decoded (> 1 => an array of values)
217 public :: fm_change_list !< (list) return success
218 public :: fm_change_root !< (list) return success
219 public :: fm_dump_list !< (list [, recursive]) return success
220 public :: fm_exists !< (field) return success
221 public :: fm_get_index !< (field) return index
222 public :: fm_get_current_list !< () return path
223 public :: fm_get_length !< (list) return length
224 public :: fm_get_type !< (field) return string
225 public :: fm_get_value !< (entry, value [, index]) return success !! generic
226 public :: fm_get_value_integer !< as above (overloaded function)
227 public :: fm_get_value_logical !< as above (overloaded function)
228 public :: fm_get_value_real_r4 !< as above (overloaded function)
229 public :: fm_get_value_real_r8 !< as above (overloaded function)
230 public :: fm_get_value_string !< as above (overloaded function)
231 public :: fm_init_loop !< (list, iter)
232 public :: fm_loop_over_list !< (list, name, type, index) return success
233 !! (iter, name, type, index) return success
234 public :: fm_new_list !< (list [, create] [, keep]) return index
235 public :: fm_new_value !< (entry, value [, create] [, index]) return index !! generic
236 public :: fm_new_value_integer !< as above (overloaded function)
237 public :: fm_new_value_logical !< as above (overloaded function)
238 public :: fm_new_value_real_r4 !< as above (overloaded function)
239 public :: fm_new_value_real_r8 !< as above (overloaded function)
240 public :: fm_new_value_string !< as above (overloaded function)
241 public :: fm_reset_loop !< ()
242 public :: fm_return_root !< () return success
243 public :: fm_modify_name !< (oldname, newname) return success
244 public :: fm_query_method !< (name, method_name, method_control) return success and
245 !! name and control strings
246 public :: fm_find_methods !< (list, methods, control) return success and name and
248 public :: fm_copy_list !< (list, suffix, [create]) return index
249 private :: create_field ! (list_p, name) return field pointer
250 private :: dump_list ! (list_p, recursive, depth) return success
251 private :: find_base ! (field, path, base)
252 private :: find_field ! (field, list_p) return field pointer
253 private :: find_head ! (field, head, rest)
254 private :: find_list ! (list, list_p, create) return field pointer
255 private :: get_field ! (field, list_p) return field pointer
256 private :: initialize_module_variables ! ()
257 private :: make_list ! (list_p, name) return field pointer
259 !> The length of a character string representing the field name.
260 integer, parameter, public :: fm_field_name_len = 48
261 !> The length of a character string representing the field path.
262 integer, parameter, public :: fm_path_name_len = 512
263 !> The length of a character string representing character values for the field.
264 integer, parameter, public :: fm_string_len = 1024
265 !> The length of a character string representing the various types that the values of the field can take.
266 integer, parameter, public :: fm_type_name_len = 8
267 !> Number of models (ATMOS, OCEAN, LAND, ICE, COUPLER).
268 integer, parameter, public :: NUM_MODELS = 5
269 !> The value returned if a field is not defined.
270 integer, parameter, public :: NO_FIELD = -1
271 !> Atmospheric model.
272 integer, parameter, public :: MODEL_ATMOS = 1
274 integer, parameter, public :: MODEL_OCEAN = 2
276 integer, parameter, public :: MODEL_LAND = 3
278 integer, parameter, public :: MODEL_ICE = 4
280 integer, parameter, public :: MODEL_COUPLER = 5
281 !> Model names, e.g. MODEL_NAMES(MODEL_OCEAN) is 'oceanic'
282 character(len=11), parameter, public, dimension(NUM_MODELS) :: &
283 MODEL_NAMES=(/'atmospheric','oceanic ','land ','ice ','coupler '/)
287 !> @brief This method_type is a way to allow a component module to alter the parameters it needs
288 !! for various tracers.
290 !> In essence this is a way to modify a namelist. A namelist can supply
291 !! default parameters for all tracers. This method will allow the user to modify these
292 !! default parameters for an individual tracer. An example could be that the user wishes to
293 !! use second order advection on a tracer and also use fourth order advection on a second
294 !! tracer within the same model run. The default advection could be second order and the
295 !! field table would then indicate that the second tracer requires fourth order advection.
296 !! This would be parsed by the advection routine.
297 !> @ingroup field_manager_mod
298 type, public :: method_type
300 character(len=fm_string_len) :: method_type !< This string represents a tag that a module
301 !! using this method can key on. Typically this should
302 !! contain some reference to the module that is calling it.
303 character(len=fm_string_len) :: method_name !< This is the name of a method which the module
304 !! can parse and use to assign different default values to
306 character(len=fm_string_len) :: method_control !< This is the string containing parameters that
307 !! the module can use as values for a field method. These should
308 !! override default values within the module.
311 !> This method_type is the same as method_type except that the
312 !! method_control string is not present. This is used when you wish to
313 !! change to a scheme within a module but do not need to pass
314 !! parameters. See @ref method_type for member information.
315 !> @ingroup field_manager_mod
316 type, public :: method_type_short
317 character(len=fm_string_len) :: method_type
318 character(len=fm_string_len) :: method_name
321 !> This is the same as method_type except that the
322 !! method_control and method_name strings are not present. This is used
323 !! when you wish to change to a scheme within a module but do not need
324 !! to pass parameters.
325 !> @ingroup field_manager_mod
326 type, public :: method_type_very_short
327 character(len=fm_string_len) :: method_type
330 !> Iterator over the field manager list
331 !> @ingroup field_manager_mod
332 type, public :: fm_list_iter_type
333 type(field_def), pointer :: ptr => NULL() !< pointer to the current field
334 end type fm_list_iter_type
336 !> @ingroup field_manager_mod
337 type(method_type), public :: default_method
339 !> @brief Returns an index corresponding to the given field name.
341 !> Model number can be given for old method.
342 !! <br>Example usage:
344 !! value=find_field_index( model, field_name )
345 !! value=find_field_index( field_name )
347 !> @ingroup field_manager_mod
348 interface find_field_index
349 module procedure find_field_index_old
350 module procedure find_field_index_new
353 !> @brief A function to parse an integer or an array of integers,
354 !! a real or an array of reals, a string or an array of strings.
356 !> Parse is an integer function that decodes values from a text string.
357 !! The text string has the form: "label=list" where "label" is an
358 !! arbitrary user defined label describing the values being decoded,
359 !! and "list" is a list of one or more values separated by commas.
360 !! The values may be integer, real, or character.
361 !! Parse returns the number of values decoded.
362 !! <br>Example usage:
364 !! number = parse(text, label, value)
366 !> @ingroup field_manager_mod
368 module procedure parse_real_r4
369 module procedure parse_real_r8
370 module procedure parse_reals_r4
371 module procedure parse_reals_r8
372 module procedure parse_integer
373 module procedure parse_integers
374 module procedure parse_string
375 module procedure parse_strings
378 !> @brief An overloaded function to assign a value to a field.
380 !> Allocate and initialize a new value and return the index.
381 !! If an error condition occurs the parameter NO_FIELD is returned.
383 !! If the type of the field is changing (e.g. real values being transformed to
384 !! integers), then any previous values for the field are removed and replaced
385 !! by the value passed in the present call to this function.
387 !! If append is present and .true., then index cannot be greater than 0 if
389 !! <br> Example usage:
391 !! field_index= fm_new_value(name, value, [create], [index], [append])
393 !> @ingroup field_manager_mod
394 interface fm_new_value
395 module procedure fm_new_value_integer
396 module procedure fm_new_value_logical
397 module procedure fm_new_value_real_r4
398 module procedure fm_new_value_real_r8
399 module procedure fm_new_value_string
402 !> @brief An overloaded function to find and extract a value for a named field.
404 !> Find and extract the value for name. The value may be of type real,
405 !! integer, logical or character. If a single value from an array of values
406 !! is required, an optional index can be supplied.
407 !! Return true for success and false for failure
408 !! <br> Example usage:
410 !! success = fm_get_value(name, value, index)
412 !> @ingroup field_manager_mod
413 interface fm_get_value
414 module procedure fm_get_value_integer
415 module procedure fm_get_value_logical
416 module procedure fm_get_value_real_r4
417 module procedure fm_get_value_real_r8
418 module procedure fm_get_value_string
421 !> @brief A function for looping over a list.
423 !> Loop over the list, setting the name, type and index
424 !! of the next field. Return false at the end of the loop.
425 !! <br> Example usage:
427 !! success = fm_loop_over_list(list, name, field_type, index)
429 !> @ingroup field_manager_mod
430 interface fm_loop_over_list
431 module procedure fm_loop_over_list_new
432 module procedure fm_loop_over_list_old
435 character(len=17), parameter :: module_name = 'field_manager_mod'
436 character(len=33), parameter :: error_header = '==>Error from '//trim(module_name)//': '
437 character(len=35), parameter :: warn_header = '==>Warning from '//trim(module_name)//': '
438 character(len=32), parameter :: note_header = '==>Note from '//trim(module_name)//': '
439 character(len=1), parameter :: comma = ","
440 character(len=1), parameter :: list_sep = '/'
441 !TODO these variable can be removed when the legacy table is no longer used
442 character(len=1), parameter :: comment = '#'
443 character(len=1), parameter :: dquote = '"'
444 character(len=1), parameter :: equal = '='
445 character(len=1), parameter :: squote = "'"
447 integer, parameter :: null_type = 0
448 integer, parameter :: integer_type = 1
449 integer, parameter :: list_type = 2
450 integer, parameter :: logical_type = 3
451 integer, parameter :: real_type = 4
452 integer, parameter :: string_type = 5
453 integer, parameter :: num_types = 5
454 integer, parameter :: array_increment = 10
455 !TODO these variable can be removed when the legacy table is no longer used
456 integer, parameter :: MAX_FIELDS = MAXFIELDS_
457 integer, parameter :: MAX_FIELD_METHODS = MAXFIELDMETHODS_
460 !> @brief Private type for internal use
461 !> @ingroup field_manager_mod
462 type, private :: field_mgr_type
463 character(len=fm_field_name_len) :: field_type
464 character(len=fm_string_len) :: field_name
465 integer :: model, num_methods
466 type(method_type), dimension(:), allocatable :: methods !< methods associated with this field name
467 end type field_mgr_type
469 !TODO These two types: field_names_type and field_names_type_short
470 !! will no longer be needed when the legacy field table is not used
471 !> @brief Private type for internal use
472 !> @ingroup field_manager_mod
473 type, private :: field_names_type
474 character(len=fm_field_name_len) :: fld_type
475 character(len=fm_field_name_len) :: mod_name
476 character(len=fm_string_len) :: fld_name
477 end type field_names_type
479 !> @brief Private type for internal use
480 !> @ingroup field_manager_mod
481 type, private :: field_names_type_short
482 character(len=fm_field_name_len) :: fld_type
483 character(len=fm_field_name_len) :: mod_name
484 end type field_names_type_short
486 !> @brief Private type for internal use
487 !> @ingroup field_manager_mod
488 type, private :: field_def
489 character (len=fm_field_name_len) :: name
491 type (field_def), pointer :: parent => NULL()
492 integer :: field_type
496 type (field_def), pointer :: first_field => NULL()
497 type (field_def), pointer :: last_field => NULL()
498 integer, allocatable, dimension(:) :: i_value
499 logical, allocatable, dimension(:) :: l_value
500 real(r8_kind), allocatable, dimension(:) :: r_value !< string to real conversion will be done at r8;
501 !! all real values will be stored as r8_kind.
502 character(len=fm_string_len), allocatable, dimension(:) :: s_value
503 type (field_def), pointer :: next => NULL()
504 type (field_def), pointer :: prev => NULL()
507 !> @addtogroup field_manager_mod
510 type(field_mgr_type), dimension(:), allocatable, private :: fields !< fields of field_mgr_type
512 character(len=fm_path_name_len) :: loop_list
513 character(len=fm_type_name_len) :: field_type_name(num_types)
514 character(len=fm_field_name_len) :: save_root_name
515 ! The string set is the set of characters.
516 character(len=52) :: set = "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz"
517 ! If a character in the string being parsed matches a character within
518 ! the string set_nonexp then the string being parsed cannot be a number.
519 character(len=50) :: set_nonexp = "ABCDFGHIJKLMNOPQRSTUVWXYZabcdfghijklmnopqrstuvwxyz"
520 ! If a character in the string being parsed matches a character within
521 ! the string setnum then the string may be a number.
522 character(len=13) :: setnum = "0123456789+-."
523 integer :: num_fields = 0
524 type (field_def), pointer :: loop_list_p => NULL()
525 type (field_def), pointer :: current_list_p => NULL()
526 type (field_def), pointer :: root_p => NULL()
527 type (field_def), pointer :: save_root_parent_p => NULL()
528 type (field_def), target, save :: root
530 logical :: use_field_table_yaml = .false. !< .True. if using the field_table.yaml,
531 !! .false. if using the legacy field_table
533 namelist /field_manager_nml/ use_field_table_yaml
537 !> @brief Routine to initialize the field manager.
539 !> This routine reads from a file containing yaml paramaters.
540 !! These yaml parameters contain information on which schemes are
541 !! needed within various modules. The field manager does not
542 !! initialize any of those schemes however. It simply holds the
543 !! information and is queried by the appropriate module.
545 !! The routine has two loops. The first loop initializes the my_table object
546 !! and counts the number of fields contained therein. The second loop is the
547 !! main loop that acts on each field in the my_table object, defining a list
548 !! object (in the field_manager definition) from which various fm routines may be
549 !! called, as well as populating the "fields" object and the "methods" objects
550 !! within each field object. The "fields" and "methods" objects are then used
551 !! with the subroutine new_name to append various characteristics to the list
552 !! object. Note that the "fields" and "methods" objects are also used with other
553 !! fm routines in a bit of a parallel system.
554 subroutine field_manager_init(nfields, table_name)
555 integer, intent(out), optional :: nfields !< number of fields
556 character(len=fm_string_len), intent(in), optional :: table_name !< Name of the field table, default
558 if (module_is_initialized) then
559 if(present(nfields)) nfields = num_fields
563 call initialize_module_variables()
565 !TODO the use_field_table_yaml namelist can be removed when the legacy table is no longer in used
566 if (use_field_table_yaml) then
567 !Crash if you are not compiling with -Duse_yaml or if the field_table is present
569 call mpp_error(FATAL, "You cannot have use_field_table_yaml=.true. without compiling with -Duse_yaml")
571 if (file_exists("field_table")) &
572 call mpp_error(FATAL, "You cannot have the legacy field_table if use_field_table_yaml=.true.")
574 call mpp_error(NOTE, "You are using the yaml version of the field_table")
575 call read_field_table_yaml(nfields, table_name)
578 if (file_exists("field_table.yaml")) &
579 call mpp_error(FATAL, "You cannot have the yaml field_table if use_field_table_yaml=.false.")
580 call mpp_error(NOTE, "You are using the legacy version of the field_table")
581 call read_field_table_legacy(nfields, table_name)
584 end subroutine field_manager_init
588 !> @brief Routine to read and parse the field table yaml
589 subroutine read_field_table_yaml(nfields, table_name)
590 integer, intent(out), optional :: nfields !< number of fields
591 character(len=fm_string_len), intent(in), optional :: table_name !< Name of the field table, default
593 character(len=fm_string_len) :: tbl_name !< field_table yaml file
594 character(len=fm_string_len) :: method_control !< field_table yaml file
595 integer :: h, i, j, k, l, m !< dummy integer buffer
596 type (fmTable_t) :: my_table !< the field table
597 integer :: model !< model assocaited with the current field
598 character(len=fm_path_name_len) :: list_name !< field_manager list name
599 character(len=fm_string_len) :: subparamvalue !< subparam value to be used when defining new name
600 character(len=fm_string_len) :: fm_yaml_null !< useful hack when OG subparam does not contain an equals sign
601 integer :: current_field !< field index within loop
602 integer :: index_list_name !< integer used as check for "no field"
603 integer :: subparamindex !< index to identify whether subparams exist for this field
604 logical :: fm_success !< logical for whether fm_change_list was a success
605 logical :: subparams !< logical whether subparams exist in this iteration
607 if (.not.PRESENT(table_name)) then
608 tbl_name = 'field_table.yaml'
610 tbl_name = trim(table_name)
612 if (.not. file_exists(trim(tbl_name))) then
613 if(present(nfields)) nfields = 0
618 ! Define my_table object and read in number of fields
619 my_table = fmTable_t(trim(tbl_name))
620 call my_table%get_blocks
621 call my_table%create_children
622 do h=1,my_table%nchildren
623 do i=1,my_table%children(h)%nchildren
624 do j=1,my_table%children(h)%children(i)%nchildren
625 num_fields = num_fields + 1
630 allocate(fields(num_fields))
633 do h=1,my_table%nchildren
634 do i=1,my_table%children(h)%nchildren
635 select case (my_table%children(h)%children(i)%name)
637 model = MODEL_COUPLER
647 call mpp_error(FATAL, trim(error_header)//'The model name is unrecognised : &
648 &'//trim(my_table%children(h)%children(i)%name))
650 do j=1,my_table%children(h)%children(i)%nchildren
651 current_field = current_field + 1
652 list_name = list_sep//lowercase(trim(my_table%children(h)%children(i)%name))//list_sep//&
653 lowercase(trim(my_table%children(h)%name))//list_sep//&
654 lowercase(trim(my_table%children(h)%children(i)%children(j)%name))
655 index_list_name = fm_new_list(list_name, create = .true.)
656 if ( index_list_name == NO_FIELD ) &
657 call mpp_error(FATAL, trim(error_header)//'Could not set field list for '//trim(list_name))
658 fm_success = fm_change_list(list_name)
659 fields(current_field)%model = model
660 fields(current_field)%field_name = lowercase(trim(my_table%children(h)%children(i)%children(j)%name))
661 fields(current_field)%field_type = lowercase(trim(my_table%children(h)%name))
662 fields(current_field)%num_methods = size(my_table%children(h)%children(i)%children(j)%key_ids)
663 allocate(fields(current_field)%methods(fields(current_field)%num_methods))
664 if(fields(current_field)%num_methods.gt.0) then
665 if (my_table%children(h)%children(i)%children(j)%nchildren .gt. 0) subparams = .true.
666 do k=1,size(my_table%children(h)%children(i)%children(j)%keys)
667 fields(current_field)%methods(k)%method_type = &
668 lowercase(trim(my_table%children(h)%children(i)%children(j)%keys(k)))
669 fields(current_field)%methods(k)%method_name = &
670 lowercase(trim(my_table%children(h)%children(i)%children(j)%values(k)))
671 if (.not.subparams) then
672 call new_name_yaml(list_name, my_table%children(h)%children(i)%children(j)%keys(k),&
673 my_table%children(h)%children(i)%children(j)%values(k) )
676 do l=1,my_table%children(h)%children(i)%children(j)%nchildren
677 if(lowercase(trim(my_table%children(h)%children(i)%children(j)%children(l)%paramname)).eq.&
678 lowercase(trim(fields(current_field)%methods(k)%method_type))) then
683 if (subparamindex.eq.-1) then
684 call new_name_yaml(list_name, my_table%children(h)%children(i)%children(j)%keys(k),&
685 my_table%children(h)%children(i)%children(j)%values(k) )
687 do m=1,size(my_table%children(h)%children(i)%children(j)%children(subparamindex)%keys)
690 if (trim(my_table%children(h)%children(i)%children(j)%values(k)).eq.'fm_yaml_null') then
693 fm_yaml_null = trim(my_table%children(h)%children(i)%children(j)%values(k))//'/'
695 method_control = trim(my_table%children(h)%children(i)%children(j)%keys(k))//"/"//&
696 &trim(fm_yaml_null)//&
697 &trim(my_table%children(h)%children(i)%children(j)%children(subparamindex)%keys(m))
698 subparamvalue = trim(my_table%children(h)%children(i)%children(j)%children(subparamindex)%values(m))
699 call new_name_yaml(list_name, method_control, subparamvalue)
709 if (present(nfields)) nfields = num_fields
710 call my_table%destruct
711 end subroutine read_field_table_yaml
713 !> @brief Subroutine to add new values to list parameters.
715 !> This subroutine uses input strings list_name, method_name
716 !! and val_name_in to add new values to the list. Given
717 !! list_name a new list item is created that is named
718 !! method_name and is given the value or values in
719 !! val_name_in. If there is more than 1 value in
720 !! val_name_in, these values should be comma-separated.
721 subroutine new_name_yaml ( list_name, method_name_in , val_name_in)
722 character(len=*), intent(in) :: list_name !< The name of the field that is of interest here.
723 character(len=*), intent(in) :: method_name_in !< The name of the method that values are
724 !! being supplied for.
725 character(len=*), intent(inout) :: val_name_in !< The value or values that will be parsed and
726 !! used as the value when creating a new field or fields.
728 character(len=fm_string_len) :: method_name !< name of method to be attached to new list
729 character(len=fm_string_len) :: val_name !< value name (to be converted to appropriate type)
730 integer, dimension(:), allocatable :: end_val !< end values in comma separated list
731 integer, dimension(:), allocatable :: start_val !< start values in comma separated list
732 integer :: i !< loop index
733 integer :: index_t !< appending index
734 integer :: num_elem !< number of elements in comma list
735 integer :: val_int !< value when converted to integer
736 integer :: val_type !< value type represented as integer for use in select case
737 logical :: append_new !< whether or not to append to existing list structure
738 logical :: val_logic !< value when converted to logical
739 real(r8_kind) :: val_real !< value when converted to real.
740 !! All strings will be converted to r8_kind reals.
742 call strip_front_blanks(val_name_in)
743 method_name = trim(method_name_in)
744 call strip_front_blanks(method_name)
750 ! If the array of values being passed in is a comma delimited list then count
751 ! the number of elements.
753 do i = 1, len_trim(val_name_in)
754 if ( val_name_in(i:i) == comma ) then
755 num_elem = num_elem + 1
759 allocate(start_val(num_elem))
760 allocate(end_val(num_elem))
762 end_val(:) = len_trim(val_name_in)
765 do i = 1, len_trim(val_name_in)
766 if ( val_name_in(i:i) == comma ) then
767 end_val(num_elem) = i-1
768 start_val(num_elem+1) = i+1
769 num_elem = num_elem + 1
775 if ( i .gt. 1 .or. index_t .eq. 0 ) then
777 index_t = 0 ! If append is true then index must be <= 0
779 val_type = string_type ! Assume it is a string
780 val_name = val_name_in(start_val(i):end_val(i))
781 call strip_front_blanks(val_name)
783 if ( scan(val_name(1:1), setnum ) > 0 ) then
784 if ( scan(val_name, set_nonexp ) .le. 0 ) then
785 if ( scan(val_name, '.') > 0 .or. scan(val_name, 'e') > 0 .or. scan(val_name, 'E') > 0) then
786 read(val_name, *) val_real
789 read(val_name, *) val_int
790 val_type = integer_type
795 if ( len_trim(val_name) == 1 .or. len_trim(val_name) == 3) then
796 if ( val_name == 't' .or. val_name == 'T' .or. val_name == '.t.' .or. val_name == '.T.' ) then
798 val_type = logical_type
800 if ( val_name == 'f' .or. val_name == 'F' .or. val_name == '.f.' .or. val_name == '.F.' ) then
802 val_type = logical_type
805 if ( trim(lowercase(val_name)) == 'true' .or. trim(lowercase(val_name)) == '.true.' ) then
807 val_type = logical_type
809 if ( trim(lowercase(val_name)) == 'false' .or. trim(lowercase(val_name)) == '.false.' ) then
811 val_type = logical_type
814 select case(val_type)
817 if ( fm_new_value( method_name, val_int, create = .true., index = index_t, append = append_new ) < 0 ) &
818 call mpp_error(FATAL, trim(error_header)//'Could not set "' // trim(val_name) // '" for '//trim(method_name)//&
819 ' (I) for '//trim(list_name))
822 if ( fm_new_value( method_name, val_logic, create = .true., index = index_t, append = append_new) < 0 ) &
823 call mpp_error(FATAL, trim(error_header)//'Could not set "' // trim(val_name) // '" for '//trim(method_name)//&
824 ' (L) for '//trim(list_name))
827 if ( fm_new_value( method_name, val_real, create = .true., index = index_t, append = append_new) < 0 ) &
828 call mpp_error(FATAL, trim(error_header)//'Could not set "' // trim(val_name) // '" for '//trim(method_name)//&
829 ' (R) for '//trim(list_name))
832 if ( fm_new_value( method_name, val_name, create = .true., index = index_t, append = append_new) < 0 ) &
833 call mpp_error(FATAL, trim(error_header)//'Could not set "' // trim(val_name) // '" for '//trim(method_name)//&
834 ' (S) for '//trim(list_name))
836 call mpp_error(FATAL, trim(error_header)//'Could not find a valid type to set the '//trim(method_name)//&
837 ' for '//trim(list_name))
842 deallocate(start_val)
845 end subroutine new_name_yaml
848 !> @brief Routine to read and parse the field table yaml
850 !> This routine reads from a file containing formatted strings.
851 !! These formatted strings contain information on which schemes are
852 !! needed within various modules. The field manager does not
853 !! initialize any of those schemes however. It simply holds the
854 !! information and is queried by the appropriate module.
855 subroutine read_field_table_legacy(nfields, table_name)
857 integer, intent(out), optional :: nfields !< number of fields
858 character(len=fm_string_len), intent(in), optional :: table_name !< Name of the field table, default
861 character(len=1024) :: record
862 character(len=fm_string_len) :: control_str
863 character(len=fm_path_name_len) :: list_name
864 character(len=fm_string_len) :: method_name
865 character(len=fm_string_len) :: name_str
866 character(len=fm_string_len) :: type_str
867 character(len=fm_string_len) :: val_name
868 character(len=fm_string_len) :: tbl_name
869 integer :: control_array(MAX_FIELDS,3)
872 integer :: index_list_name
882 logical :: flag_method
883 logical :: fm_success
884 type(field_names_type_short) :: text_names_short
885 type(field_names_type) :: text_names
886 type(method_type_short) :: text_method_short
887 type(method_type) :: text_method
888 type(method_type_very_short) :: text_method_very_short
890 if (.not.PRESENT(table_name)) then
891 tbl_name = 'field_table'
893 tbl_name = trim(table_name)
895 if (.not. file_exists(trim(tbl_name))) then
896 if(present(nfields)) nfields = 0
900 allocate(fields(MAX_FIELDS))
902 open(newunit=iunit, file=trim(tbl_name), action='READ', iostat=io_status)
903 if(io_status/=0) call mpp_error(FATAL, 'field_manager_mod: Error in opening file '//trim(tbl_name))
904 !write_version_number should precede all writes to stdlog from field_manager
905 call write_version_number("FIELD_MANAGER_MOD", version)
908 read(iunit,'(a)',end=89,err=99) record
909 write( log_unit,'(a)' )record
910 if (record(1:1) == "#" ) cycle
911 ltrec = LEN_TRIM(record)
912 if (ltrec .le. 0 ) cycle ! Blank line
917 if (record(l:l) == '"' ) then
921 if (icount > 6 ) then
922 call mpp_error(FATAL,trim(error_header)//'Too many fields in field table header entry.'//trim(record))
927 read(record,*,end=79,err=79) text_names
928 text_names%fld_type = lowercase(trim(text_names%fld_type))
929 text_names%mod_name = lowercase(trim(text_names%mod_name))
930 text_names%fld_name = lowercase(trim(text_names%fld_name))
932 ! If there is no control string then the last string can be omitted and there are only 4 '"' in the record.
933 read(record,*,end=79,err=79) text_names_short
934 text_names%fld_type = lowercase(trim(text_names_short%fld_type))
935 text_names%mod_name = lowercase(trim(text_names_short%mod_name))
936 text_names%fld_name = lowercase(trim(text_names_short%mod_name))
938 ! If there is only the method_type string then the last 2 strings need to be blank and there
939 ! are only 2 '"' in the record.
940 read(record,*,end=79,err=79) text_names_short
941 text_names%fld_type = lowercase(trim(text_names_short%fld_type))
942 text_names%mod_name = lowercase(trim(text_names_short%mod_name))
943 text_names%fld_name = lowercase(trim(text_names_short%mod_name))
945 ! There is an unterminated or unquoted string in the field table entry.
946 text_names%fld_type = " "
947 text_names%mod_name = lowercase(trim(record))
948 text_names%fld_name = " "
951 ! Create a list with Rick Slaters field manager code
953 list_name = list_sep//trim(text_names%mod_name)//list_sep//trim(text_names%fld_type)//&
954 list_sep//trim(text_names%fld_name)
956 index_list_name = fm_new_list(list_name, create = .true.)
957 if ( index_list_name == NO_FIELD ) &
958 call mpp_error(FATAL, trim(error_header)//'Could not set field list for '//trim(list_name))
960 fm_success = fm_change_list(list_name)
961 select case (text_names%mod_name)
963 model = MODEL_COUPLER
973 call mpp_error(FATAL, trim(error_header)//'The model name is unrecognised : '//trim(text_names%mod_name))
975 if (find_field_index(list_name) > 0) then
976 num_fields = num_fields + 1
978 if (num_fields > MAX_FIELDS) call mpp_error(FATAL,trim(error_header)//'max fields exceeded')
979 fields(num_fields)%model = model
980 fields(num_fields)%field_name = lowercase(trim(text_names%fld_name))
981 fields(num_fields)%field_type = lowercase(trim(text_names%fld_type))
982 fields(num_fields)%num_methods = 0
983 allocate(fields(num_fields)%methods(MAX_FIELD_METHODS))
984 call check_for_name_duplication
986 ! Check to see that the first line is not the only line
987 if ( record(LEN_TRIM(record):LEN_TRIM(record)) == list_sep) cycle
991 do while (flag_method)
992 read(iunit,'(a)',end=99,err=99) record
993 ! If the line is blank then fetch the next line.
994 if (LEN_TRIM(record) .le. 0) cycle
995 ! If the last character in the line is / then this is the end of the field methods
996 if ( record(LEN_TRIM(record):LEN_TRIM(record)) == list_sep) then
997 flag_method = .FALSE.
998 if (LEN_TRIM(record) == 1) cycle
999 record = record(:LEN_TRIM(record)-1) ! Remove the end of field method marker
1001 ! If the line is now blank, after removing the field separator marker, then fetch the next line.
1002 if (LEN_TRIM(record) .le. 0) cycle
1003 ! If the first character in the line is # then it is treated as a comment
1004 if (record(1:1) == comment ) cycle
1007 do l= 1, LEN_TRIM(record)
1008 if (record(l:l) == dquote ) then
1012 if (icount > 6 ) call mpp_error(FATAL,trim(error_header)//'Too many fields in field entry.'//trim(record))
1014 if (.not. fm_change_list ( list_name)) &
1015 call mpp_error(FATAL, trim(error_header)//'Could not change to '//trim(list_name)//' list')
1017 select case (icount)
1019 read(record,*,end=99,err=99) text_method
1020 fields(num_fields)%methods(m)%method_type = lowercase(trim(text_method%method_type))
1021 fields(num_fields)%methods(m)%method_name = lowercase(trim(text_method%method_name))
1022 fields(num_fields)%methods(m)%method_control = lowercase(trim(text_method%method_control))
1024 type_str = text_method%method_type
1025 name_str = text_method%method_name
1026 control_str = text_method%method_control
1029 ! If there is no control string then the last string can be omitted and there are only 4 '"' in the record.
1030 read(record,*,end=99,err=99) text_method_short
1031 fields(num_fields)%methods(m)%method_type =&
1032 & lowercase(trim(text_method_short%method_type))
1033 fields(num_fields)%methods(m)%method_name =&
1034 & lowercase(trim(text_method_short%method_name))
1035 fields(num_fields)%methods(m)%method_control = " "
1037 type_str = text_method_short%method_type
1039 control_str = text_method_short%method_name
1042 ! If there is only the method_type string then the last 2 strings need to be blank and there
1043 ! are only 2 '"' in the record.
1044 read(record,*,end=99,err=99) text_method_very_short
1045 fields(num_fields)%methods(m)%method_type = lowercase(trim(text_method_very_short%method_type))
1046 fields(num_fields)%methods(m)%method_name = " "
1047 fields(num_fields)%methods(m)%method_control = " "
1051 control_str = text_method_very_short%method_type
1054 read(record,'(A)',end=99,err=99) control_str
1059 call mpp_error(FATAL,trim(error_header)//'Unterminated field in field entry.'//trim(record))
1062 ! This section of code breaks the control string into separate strings.
1063 ! The array control_array contains the following parameters.
1064 ! control_array(:,1) = index within control_str of the first character of the name.
1065 ! control_array(:,2) = index within control_str of the equal sign
1066 ! control_array(:,3) = index within control_str of the last character of the value.
1068 ! control_array(:,1) -> control_array(:,2) -1 = name of the parameter.
1069 ! control_array(:,2)+1 -> control_array(:,3) = value of the parameter.
1071 ltrec= len_trim(control_str)
1072 control_array(:,1) = 1
1073 control_array(:,2:3) = ltrec
1076 if (control_str(l:l) == equal ) then
1078 control_array(icount,2) = l ! Middle of string
1079 elseif (control_str(l:l) == comma ) then
1080 if (icount .eq. 0) then
1081 call mpp_error(FATAL,trim(error_header) // &
1082 ' Bad format for field entry (comma without equals sign): ''' // &
1083 trim(control_str) // '''')
1084 elseif (icount .gt. MAX_FIELDS) then
1085 call mpp_error(FATAL,trim(error_header) // &
1086 ' Too many fields in field entry: ''' // &
1087 trim(control_str) // '''')
1089 control_array(icount,3) = l-1 !End of previous string
1090 control_array(min(MAX_FIELDS,icount+1),1) = l+1 !Start of next string
1095 ! Make sure that we point to the end of the string (minus any trailing comma)
1096 ! for the last set of values. This fixes the case where the last set of values
1097 ! is a comma separated list
1099 if (control_str(ltrec:ltrec) .ne. comma) then
1100 control_array(max(1,icount),3) = ltrec
1103 if ( icount == 0 ) then
1104 method_name = type_str
1105 if (len_trim(method_name) > 0 ) then
1106 method_name = trim(method_name)//list_sep// trim(name_str)
1108 method_name = trim(name_str)
1110 val_name = control_str
1112 call new_name(list_name, method_name, val_name )
1117 startcont = control_array(l,1)
1118 midcont = control_array(l,2)
1119 endcont = control_array(l,3)
1121 method_name = trim(type_str)
1122 if (len_trim(method_name) > 0 ) then
1123 method_name = trim(method_name)//list_sep// trim(name_str)
1125 method_name = trim(name_str)
1128 if (len_trim(method_name) > 0 ) then
1129 method_name = trim(method_name)//list_sep//&
1130 trim(control_str(startcont:midcont-1))
1132 method_name = trim(control_str(startcont:midcont-1))
1134 val_name = trim(control_str(midcont+1:endcont))
1136 call new_name(list_name, method_name, val_name )
1141 fields(num_fields)%num_methods = fields(num_fields)%num_methods + 1
1142 if (fields(num_fields)%num_methods > MAX_FIELD_METHODS) &
1143 call mpp_error(FATAL,trim(error_header)//'Maximum number of methods for field exceeded')
1147 flag_method = .TRUE.
1148 do while (flag_method)
1149 read(iunit,'(A)',end=99,err=99) record
1150 if ( record(LEN_TRIM(record):LEN_TRIM(record)) == list_sep) then
1151 flag_method = .FALSE.
1159 close(iunit, iostat=io_status)
1160 if(io_status/=0) call mpp_error(FATAL, 'field_manager_mod: Error in closing file '//trim(tbl_name))
1163 if(present(nfields)) nfields = num_fields
1165 default_method%method_type = 'none'
1166 default_method%method_name = 'none'
1167 default_method%method_control = 'none'
1172 call mpp_error(FATAL,trim(error_header)//' Error reading field table. Record = '//trim(record))
1174 end subroutine read_field_table_legacy
1176 subroutine check_for_name_duplication
1179 ! Check that name is unique amoung fields of the same field_type and model.
1181 if ( fields(i)%field_type == fields(num_fields)%field_type .and. &
1182 fields(i)%model == fields(num_fields)%model .and. &
1183 fields(i)%field_name == fields(num_fields)%field_name ) then
1184 if (mpp_pe() .eq. mpp_root_pe()) then
1185 call mpp_error(WARNING,'Error in field_manager_mod. Duplicate field name: Field type='//&
1186 trim(fields(i)%field_type)// &
1187 ', Model='//trim(MODEL_NAMES(fields(i)%model))// &
1188 ', Duplicated name='//trim(fields(i)%field_name))
1193 end subroutine check_for_name_duplication
1195 !> @brief Subroutine to add new values to list parameters.
1197 !> This subroutine uses input strings list_name, method_name
1198 !! and val_name_in to add new values to the list. Given
1199 !! list_name a new list item is created that is named
1200 !! method_name and is given the value or values in
1201 !! val_name_in. If there is more than 1 value in
1202 !! val_name_in, these values should be comma-separated.
1203 subroutine new_name ( list_name, method_name_in , val_name_in)
1204 character(len=*), intent(in) :: list_name !< The name of the field that is of interest here.
1205 character(len=*), intent(in) :: method_name_in !< The name of the method that values are
1206 !! being supplied for.
1207 character(len=*), intent(inout) :: val_name_in !< The value or values that will be parsed and
1208 !! used as the value when creating a new field or fields.
1210 character(len=fm_string_len) :: method_name
1211 character(len=fm_string_len) :: val_name
1212 integer, dimension(MAX_FIELDS) :: end_val
1213 integer, dimension(MAX_FIELDS) :: start_val
1222 logical :: append_new
1223 logical :: val_logic
1224 real(r8_kind) :: val_real !< all reals converted from string will be in r8_kind precision
1227 call strip_front_blanks(val_name_in)
1228 method_name = trim (method_name_in)
1229 call strip_front_blanks(method_name)
1233 append_new = .false.
1235 end_val(:) = len_trim(val_name_in)
1237 ! If the array of values being passed in is a comma delimited list then count
1238 ! the number of elements.
1240 do i = 1, len_trim(val_name_in)
1241 if ( val_name_in(i:i) == comma ) then
1242 end_val(num_elem) = i-1
1243 start_val(num_elem+1) = i+1
1244 num_elem = num_elem + 1
1248 ! Check to see if this is an array element of form array[x] = value
1249 left_br = scan(method_name,'[')
1250 right_br = scan(method_name,']')
1251 if ( num_elem .eq. 1 ) then
1252 if ( left_br > 0 .and. right_br == 0 ) &
1253 call mpp_error(FATAL, trim(error_header)//"Left bracket present without right bracket in "//trim(method_name))
1254 if ( left_br== 0 .and. right_br > 0 ) &
1255 call mpp_error(FATAL, trim(error_header)//"Right bracket present without left bracket in "//trim(method_name))
1256 if ( left_br > 0 .and. right_br > 0 ) then
1257 if ( scan( method_name(left_br+1:right_br -1), set ) > 0 ) &
1258 call mpp_error(FATAL, trim(error_header)//"Using a non-numeric value for index in "//trim(method_name))
1259 read(method_name(left_br+1:right_br -1), *) index_t
1260 method_name = method_name(:left_br -1)
1263 ! If there are multiple values then there cannot be a bracket in method_name.
1264 if ( left_br > 0 .or. right_br > 0 ) &
1265 call mpp_error(FATAL, &
1266 trim(error_header)//"Using a comma delimited list with an indexed array element in "//trim(method_name))
1271 if ( i .gt. 1 .or. index_t .eq. 0 ) then
1273 index_t = 0 ! If append is true then index must be <= 0
1275 val_type = string_type ! Assume it is a string
1276 val_name = val_name_in(start_val(i):end_val(i))
1277 call strip_front_blanks(val_name)
1280 ! if the string starts and ends with matching single quotes, then this is a string
1281 ! if there are quotes which do not match, then this is an error
1284 length = len_trim(val_name)
1285 if (val_name(1:1) .eq. squote) then
1287 if (val_name(length:length) .eq. squote) then
1288 val_name = val_name(2:length-1)//repeat(" ",len(val_name)-length+2)
1289 val_type = string_type
1290 elseif (val_name(length:length) .eq. dquote) then
1291 call mpp_error(FATAL, trim(error_header) // ' Quotes do not match in ' // trim(val_name) // &
1292 ' for ' // trim(method_name) // ' of ' // trim(list_name))
1294 call mpp_error(FATAL, trim(error_header) // ' No trailing quote in ' // trim(val_name) // &
1295 ' for ' // trim(method_name) // ' of ' // trim(list_name))
1298 elseif (val_name(1:1) .eq. dquote .or. val_name(length:length) .eq. dquote) then
1300 call mpp_error(FATAL, trim(error_header) // ' Double quotes not allowed in ' // trim(val_name) // &
1301 ' for ' // trim(method_name) // ' of ' // trim(list_name))
1303 elseif (val_name(length:length) .eq. squote) then
1305 call mpp_error(FATAL, trim(error_header) // ' No leading quote in ' // trim(val_name) // &
1306 ' for ' // trim(method_name) // ' of ' // trim(list_name))
1309 ! If the string to be parsed is a real then all the characters must be numeric,
1310 ! be a plus/minus, be a decimal point or, for exponentials, be e or E.
1312 ! If a string is an integer, then all the characters must be numeric.
1314 if ( scan(val_name(1:1), setnum ) > 0 ) then
1316 ! If there is a letter in the name it may only be e or E
1318 if ( scan(val_name, set_nonexp ) .le. 0 ) then
1319 ! It is real if there is a . in the name or the value appears exponential
1320 if ( scan(val_name, '.') > 0 .or. scan(val_name, 'e') > 0 .or. scan(val_name, 'E') > 0) then
1321 read(val_name, *) val_real
1322 val_type = real_type
1324 read(val_name, *) val_int
1325 val_type = integer_type
1331 ! If val_name is t/T or f/F then this is a logical flag.
1332 if ( len_trim(val_name) == 1 .or. len_trim(val_name) == 3) then
1333 if ( val_name == 't' .or. val_name == 'T' .or. val_name == '.t.' .or. val_name == '.T.' ) then
1335 val_type = logical_type
1337 if ( val_name == 'f' .or. val_name == 'F' .or. val_name == '.f.' .or. val_name == '.F.' ) then
1339 val_type = logical_type
1342 if ( trim(lowercase(val_name)) == 'true' .or. trim(lowercase(val_name)) == '.true.' ) then
1344 val_type = logical_type
1346 if ( trim(lowercase(val_name)) == 'false' .or. trim(lowercase(val_name)) == '.false.' ) then
1348 val_type = logical_type
1352 select case(val_type)
1355 if ( fm_new_value( method_name, val_int, create = .true., index = index_t, append = append_new ) < 0 ) &
1356 call mpp_error(FATAL, trim(error_header)//'Could not set "' // trim(val_name) // '" for '//trim(method_name)//&
1357 ' (I) for '//trim(list_name))
1360 if ( fm_new_value( method_name, val_logic, create = .true., index = index_t, append = append_new) < 0 ) &
1361 call mpp_error(FATAL, trim(error_header)//'Could not set "' // trim(val_name) // '" for '//trim(method_name)//&
1362 ' (L) for '//trim(list_name))
1365 if ( fm_new_value( method_name, val_real, create = .true., index = index_t, append = append_new) < 0 ) &
1366 call mpp_error(FATAL, trim(error_header)//'Could not set "' // trim(val_name) // '" for '//trim(method_name)//&
1367 ' (R) for '//trim(list_name))
1370 if ( fm_new_value( method_name, val_name, create = .true., index = index_t, append = append_new) < 0 ) &
1371 call mpp_error(FATAL, trim(error_header)//'Could not set "' // trim(val_name) // '" for '//trim(method_name)//&
1372 ' (S) for '//trim(list_name))
1374 call mpp_error(FATAL, trim(error_header)//'Could not find a valid type to set the '//trim(method_name)//&
1375 ' for '//trim(list_name))
1381 end subroutine new_name
1383 !> @brief Destructor for field manager.
1385 !> This subroutine deallocates allocated variables (if allocated) and
1386 !! changes the initialized flag to false.
1387 subroutine field_manager_end
1390 module_is_initialized = .false.
1393 if(allocated(fields(j)%methods)) deallocate(fields(j)%methods)
1395 if(allocated(fields)) deallocate(fields)
1397 end subroutine field_manager_end
1399 !> @brief A routine to strip whitespace from the start of character strings.
1401 !> This subroutine removes spaces and tabs from the start of a character string.
1402 subroutine strip_front_blanks(name)
1404 character(len=*), intent(inout) :: name !< name to remove whitespace from
1406 name = trim(adjustl(name))
1407 end subroutine strip_front_blanks
1409 !> @brief Function to return the index of the field
1411 !> This function when passed a model number and a field name will
1412 !! return the index of the field within the field manager. This index
1413 !! can be used to access other information from the field manager.
1414 !! @returns The index of the field corresponding to field_name.
1415 function find_field_index_old(model, field_name)
1417 integer :: find_field_index_old
1418 integer, intent(in) :: model !< The number indicating which model is used.
1419 character(len=*), intent(in) :: field_name !< The name of the field that an index is being requested for.
1423 find_field_index_old = NO_FIELD
1426 if (fields(i)%model == model .and. fields(i)%field_name == lowercase(field_name)) then
1427 find_field_index_old = i
1432 end function find_field_index_old
1434 !> @returns index of the field corresponding to field_name
1435 function find_field_index_new(field_name)
1437 integer :: find_field_index_new
1438 character(len=*), intent(in) :: field_name !< The path to the name of the field that an index is
1439 !! being requested for.
1441 find_field_index_new = NO_FIELD
1443 find_field_index_new = fm_get_index(field_name)
1445 end function find_field_index_new
1447 !> @brief This routine allows access to field information given an index.
1449 !> When passed an index, this routine will return the type of field,
1450 !! the name of the field, the model which the field is associated and
1451 !! the number of methods associated with the field.
1452 !! <br>Example usage:
1454 !! call get_field_info( n,fld_type,fld_name,model,num_methods )
1456 subroutine get_field_info(n,fld_type,fld_name,model,num_methods)
1457 integer, intent(in) :: n !< index of field
1458 character (len=*),intent(out) :: fld_type !< field type
1459 character (len=*),intent(out) :: fld_name !< name of the field
1460 integer, intent(out) :: model !< number indicating which model is used
1461 integer, intent(out) :: num_methods !< number of methods
1463 if (n < 1 .or. n > num_fields) call mpp_error(FATAL,trim(error_header)//'Invalid field index')
1465 fld_type = fields(n)%field_type
1466 fld_name = fields(n)%field_name
1467 model = fields(n)%model
1468 num_methods = fields(n)%num_methods
1470 end subroutine get_field_info
1472 !> @brief A routine to get a specified method
1474 !> This routine, when passed a field index and a method index will
1475 !! return the method text associated with the field(n) method(m).
1476 subroutine get_field_method(n,m,method)
1478 integer, intent(in) :: n !< index of field
1479 integer, intent(in) :: m !< index of method
1480 type(method_type) ,intent(inout) :: method !< the m-th method of field with index n
1482 if (n < 1 .or. n > num_fields) call mpp_error(FATAL,trim(error_header)//'Invalid field index')
1483 if (m < 1 .or. m > fields(n)%num_methods) call mpp_error(FATAL,trim(error_header)//'Invalid method index')
1485 method = fields(n)%methods(m)
1487 end subroutine get_field_method
1489 !> @brief A routine to obtain all the methods associated with a field.
1491 !> When passed a field index, this routine will return the text
1492 !! associated with all the methods attached to the field.
1493 subroutine get_field_methods(n,methods)
1495 integer, intent(in) :: n !< field index
1496 type(method_type),intent(inout) :: methods(:) !< an array of methods for field with index n
1498 if (n < 1 .or. n > num_fields) &
1499 call mpp_error(FATAL,trim(error_header)//'Invalid field index')
1501 if (size(methods(:)) < fields(n)%num_methods) &
1502 call mpp_error(FATAL,trim(error_header)//'Method array too small')
1504 methods = default_method
1505 methods(1:fields(n)%num_methods) = fields(n)%methods(1:fields(n)%num_methods)
1507 end subroutine get_field_methods
1509 !> @returns The number of values that have been decoded. This allows
1510 !! a user to define a large array and fill it partially with
1511 !! values from a list. This should be the size of the value array.
1512 function parse_integers ( text, label, values ) result (parse)
1513 character(len=*), intent(in) :: text !< The text string from which the values will be parsed.
1514 character(len=*), intent(in) :: label !< A label which describes the values being decoded.
1515 integer, intent(out) :: values(:) !< The value or values that have been decoded.
1518 end function parse_integers
1520 function parse_strings ( text, label, values ) result (parse)
1521 character(len=*), intent(in) :: text !< The text string from which the values will be parsed.
1522 character(len=*), intent(in) :: label !< A label which describes the values being decoded.
1523 character(len=*), intent(out) :: values(:) !< The value or values that have been decoded.
1526 end function parse_strings
1528 function parse_integer ( text, label, parse_ival ) result (parse)
1529 character(len=*), intent(in) :: text !< The text string from which the values will be parsed.
1530 character(len=*), intent(in) :: label !< A label which describes the values being decoded.
1531 integer, intent(out) :: parse_ival !< The value or values that have been decoded.
1534 integer :: values(1)
1536 parse = parse_integers ( text, label, values )
1537 if (parse > 0) parse_ival = values(1)
1538 end function parse_integer
1540 function parse_string ( text, label, parse_sval ) result (parse)
1541 character(len=*), intent(in) :: text !< The text string from which the values will be parsed.
1542 character(len=*), intent(in) :: label !< A label which describes the values being decoded.
1543 character(len=*), intent(out) :: parse_sval !< The value or values that have been decoded.
1546 character(len=len(parse_sval)) :: values(1)
1548 parse = parse_strings ( text, label, values )
1549 if (parse > 0) parse_sval = values(1)
1550 end function parse_string
1552 !> @brief A function to create a field as a child of parent_p. This will return
1553 !! a pointer to a field_def type.
1555 !> Allocate and initialize a new field in parent_p list.
1556 !! Return a pointer to the field on success, or a null pointer
1558 !! <br>Example usage:
1560 !! list_p => create_field(parent_p, name)
1562 function create_field(parent_p, name) &
1564 type (field_def), pointer :: list_p
1565 type (field_def), pointer :: parent_p !< A pointer to the parent of the field that is to be created
1566 character(len=*), intent(in) :: name !< The name of the field that is to be created
1568 integer :: error, out_unit
1569 ! Check for fatal errors which should never arise
1571 if (.not. associated(parent_p) .or. name .eq. ' ') then
1576 ! Allocate space for the new list
1577 allocate(list_p, stat = error)
1578 if (error .ne. 0) then
1579 write (out_unit,*) trim(error_header), 'Error ', error, &
1580 ' allocating memory for list ', trim(name)
1584 ! Initialize the new field
1587 nullify(list_p%next)
1588 list_p%prev => parent_p%last_field
1589 nullify(list_p%first_field)
1590 nullify(list_p%last_field)
1592 list_p%field_type = null_type
1593 list_p%max_index = 0
1594 list_p%array_dim = 0
1595 if (allocated(list_p%i_value)) deallocate(list_p%i_value)
1596 if (allocated(list_p%l_value)) deallocate(list_p%l_value)
1597 if (allocated(list_p%r_value)) deallocate(list_p%r_value)
1598 if (allocated(list_p%s_value)) deallocate(list_p%s_value)
1599 ! If this is the first field in the parent, then set the pointer
1600 ! to it, otherwise, update the "next" pointer for the last list
1601 if (parent_p%length .le. 0) then
1602 parent_p%first_field => list_p
1604 parent_p%last_field%next => list_p
1606 ! Update the pointer for the last list in the parent
1607 parent_p%last_field => list_p
1608 ! Update the length for the parent
1609 parent_p%length = parent_p%length + 1
1610 ! Set the new index as the return value
1611 list_p%index = parent_p%length
1612 ! set the pointer to the parent list
1613 list_p%parent => parent_p
1615 end function create_field
1617 !> @brief This is a function that lists the parameters of a field.
1619 !> Given a pointer to a list, this function prints out the fields, and
1620 !! subfields, if recursive is true, associated with the list.
1622 !! This is most likely to be used through fm_dump_list.
1623 !! <br> Example usage:
1625 !! success = dump_list(list_p, recursive=.true., depth=0)
1627 logical recursive function dump_list(list_p, recursive, depth, out_unit) result(success)
1629 type (field_def), pointer :: list_p !< pointer to the field to be printed out
1630 logical, intent(in) :: recursive !< flag to make function recursively print subfields
1631 integer, intent(in) :: depth !< Listing will be padded so that 'depth' spaces appear before
1632 !! the field being printed
1633 integer, intent(in) :: out_unit !< unit number to print to
1637 character(len=fm_field_name_len) :: num, scratch
1638 type (field_def), pointer :: this_field_p
1639 character(len=depth+fm_field_name_len) :: blank
1641 blank = ' ' ! initialize blank string
1643 ! Check for a valid list
1645 if (.not. associated(list_p)) then
1647 elseif (list_p%field_type .ne. list_type) then
1651 ! set the default return value
1654 ! Print the name of this list
1655 write (out_unit,'(a,a,a)') blank(1:depth), trim(list_p%name), list_sep
1657 ! Increment the indentation depth
1658 ! The following max function is to work around an error in the IBM compiler for len_trim
1659 ! depthp1 = depth + max(len_trim(list_p%name),0) + len_trim(list_sep)
1662 this_field_p => list_p%first_field
1664 do while (associated(this_field_p))
1666 select case(this_field_p%field_type)
1668 ! If this is a list, then call dump_list
1670 ! If recursive is true, then this routine will find and dump sub-fields.
1671 success = dump_list(this_field_p, .true., depthp1, out_unit)
1672 if (.not.success) exit ! quit immediately in case of error
1673 else ! Otherwise it will print out the name of this field.
1674 write (out_unit,'(a,a,a)') blank(1:depthp1), trim(this_field_p%name), list_sep
1678 if (this_field_p%max_index .eq. 0) then
1679 write (out_unit,'(a,a,a)') blank(1:depthp1), trim(this_field_p%name), ' = NULL'
1680 elseif (this_field_p%max_index .eq. 1) then
1681 write (scratch,*) this_field_p%i_value(1)
1682 write (out_unit,'(a,a,a,a)') blank(1:depthp1), trim(this_field_p%name), ' = ', &
1683 trim(adjustl(scratch))
1684 else ! Write out the array of values for this field.
1685 do j = 1, this_field_p%max_index
1686 write (scratch,*) this_field_p%i_value(j)
1688 write (out_unit,'(a,a,a,a,a,a)') blank(1:depthp1), trim(this_field_p%name), &
1689 '[', trim(adjustl(num)), '] = ', trim(adjustl(scratch))
1694 if (this_field_p%max_index .eq. 0) then
1695 write (out_unit,'(a,a,a)') blank(1:depthp1), trim(this_field_p%name), ' = NULL'
1696 elseif (this_field_p%max_index .eq. 1) then
1697 write (scratch,'(l1)') this_field_p%l_value(1)
1698 write (out_unit,'(a,a,a,a)') blank(1:depthp1), trim(this_field_p%name), ' = ', &
1699 trim(adjustl(scratch))
1700 else ! Write out the array of values for this field.
1701 do j = 1, this_field_p%max_index
1702 write (scratch,'(l1)') this_field_p%l_value(j)
1704 write (out_unit,'(a,a,a,a,a,a)') blank(1:depthp1), trim(this_field_p%name), &
1705 '[', trim(adjustl(num)), '] = ', trim(adjustl(scratch))
1710 if (this_field_p%max_index .eq. 0) then
1711 write (out_unit,'(a,a,a)') blank(1:depthp1), trim(this_field_p%name), ' = NULL'
1712 elseif (this_field_p%max_index .eq. 1) then
1713 write (scratch,*) this_field_p%r_value(1)
1714 write (out_unit,'(a,a,a,a)') blank(1:depthp1), trim(this_field_p%name), ' = ', &
1715 trim(adjustl(scratch))
1716 else ! Write out the array of values for this field.
1717 do j = 1, this_field_p%max_index
1718 write (scratch,*) this_field_p%r_value(j)
1720 write (out_unit,'(a,a,a,a,a,a)') blank(1:depthp1), trim(this_field_p%name), &
1721 '[', trim(adjustl(num)), '] = ', trim(adjustl(scratch))
1726 if (this_field_p%max_index .eq. 0) then
1727 write (out_unit,'(a,a,a)') blank(1:depthp1), trim(this_field_p%name), ' = NULL'
1728 elseif (this_field_p%max_index .eq. 1) then
1729 write (out_unit,'(a,a,a,a)') blank(1:depthp1), trim(this_field_p%name), ' = ', &
1730 ''''//trim(this_field_p%s_value(1))//''''
1731 else ! Write out the array of values for this field.
1732 do j = 1, this_field_p%max_index
1734 write (out_unit,'(a,a,a,a,a,a)') blank(1:depthp1), trim(this_field_p%name), &
1735 '[', trim(adjustl(num)), '] = ', ''''//trim(this_field_p%s_value(j))//''''
1745 this_field_p => this_field_p%next
1748 end function dump_list
1750 !> @brief A subroutine that splits a listname into a path and a base.
1752 !> Find the base name for a list by splitting the list name into
1753 !! a path and base. The base is the last field within name, while the
1754 !! path is the preceding section of name. The base string can then be
1755 !! used to query for values associated with name.
1756 subroutine find_base(name, path, base)
1758 character(len=*), intent(in) :: name !< list name for a field
1759 character(len=*), intent(out) :: path !< path of the base field
1760 character(len=*), intent(out) :: base !< A string which can be used to query for values associated with name
1765 ! Check for the last occurrence of the list separator in name
1766 ! The following max function is to work around an error in the IBM compiler for len_trim
1767 length = max(len_trim(name),0)
1768 if (length .eq. 0) then
1770 ! Empty name, so return empty path and base
1774 ! Remove trailing list separators
1775 do while (name(length:length) .eq. list_sep)
1777 if (length .eq. 0) then
1781 if (length .eq. 0) then
1783 ! Name only list separators, so return empty path and base
1787 ! Check for the last occurrence of the list separator in name
1788 i = index(name(1:length), list_sep, back = .true.)
1790 ! no list separators in the path, so return an empty path
1791 ! and name as the base
1793 base = name(1:length)
1795 ! Found a list separator, so return the part up to the last
1796 ! list separator in path, and the remainder in base
1798 base = name(i+1:length)
1803 end subroutine find_base
1805 !> @brief Find and return a pointer to the field in the specified
1806 !! list. Return a null pointer on error.
1808 !> Find and return a pointer to the field in the specified
1809 !! list. Return a null pointer on error. Given a pointer to a field,
1810 !! this function searchs for "name" as a sub field.
1811 !> @returns A pointer to the field corresponding to "name" or an unassociated pointer if the field
1812 !! name does not exist.
1813 function find_field(name, this_list_p) &
1815 type (field_def), pointer :: field_p
1816 character(len=*), intent(in) :: name !< The name of a field that the user wishes to find
1817 type (field_def), pointer :: this_list_p !< A pointer to a list which the user wishes to search
1818 !! for a field "name".
1820 type (field_def), pointer, save :: temp_p
1825 if (name .eq. '.') then
1827 ! If the field is '.' then return this list
1828 field_p => this_list_p
1829 elseif (name .eq. '..') then
1830 ! If the field is '..' then return the parent list
1831 field_p => this_list_p%parent
1833 ! Loop over each field in this list
1834 temp_p => this_list_p%first_field
1836 do while (associated(temp_p))
1837 ! If the name matches, then set the return pointer and exit
1839 if (temp_p%name .eq. name) then
1844 temp_p => temp_p%next
1849 end function find_field
1851 !> @brief Find the first list for a name by splitting the name into
1852 !! a head and the rest.
1854 !> Find the first list for a name by splitting the name into a head and the
1855 !! rest. The head is the first field within name, while rest is the remaining
1856 !! section of name. The head string can then be used to find other fields that
1857 !! may be associated with name.
1858 subroutine find_head(name, head, rest)
1860 character(len=*), intent(in) :: name !< The name of a field of interest
1861 character(len=*), intent(out) :: head !< the first field within name
1862 character(len=*), intent(out) :: rest !< the remaining section of name
1865 ! Check for the first occurrence of the list separator in name
1866 i = index(name, list_sep)
1867 ! Check for additional consecutive list separators and return
1869 do while (i .le. len(name))
1870 if (name(i+1:i+1) .eq. list_sep) then
1878 ! no list separators in the path, so return an empty head and
1882 elseif (i .eq. len(name)) then
1883 ! The last character in name is a list separator, so return name
1884 ! as head and an empty rest
1888 ! Found a list separator, so return the part up to the list
1889 ! separator in head, and the remainder in rest
1894 end subroutine find_head
1896 !> @brief Find and return a pointer to the specified list, relative to
1897 !! relative_p. Return a null pointer on error.
1899 !> This function, when supplied a pointer to a field and a name of a second
1900 !! field relative to that pointer, will find a list and return the pointer to
1901 !! the second field. If create is .true. and the second field does not exist,
1902 !! it will be created.
1903 !> @returns A pointer to the list to be returned
1904 function find_list(path, relative_p, create) &
1906 type (field_def), pointer :: list_p
1907 character(len=*), intent(in) :: path !< path to the list of interest
1908 type (field_def), pointer :: relative_p !< pointer to the list to which "path" is relative to
1909 logical, intent(in) :: create !< If the list does not exist, it will be created if set to true
1911 character(len=fm_path_name_len) :: working_path
1912 character(len=fm_path_name_len) :: rest
1913 character(len=fm_field_name_len) :: this_list
1914 integer :: i, out_unit
1915 type (field_def), pointer, save :: working_path_p
1916 type (field_def), pointer, save :: this_list_p
1920 ! If the path is empty, then return the relative list
1921 if (path .eq. ' ') then
1923 list_p => relative_p
1926 ! If a fully qualified path is given (i.e., starts with the
1927 ! list separator) then do everything relative to root,
1928 ! otherwise, do everything relative to relative list.
1929 if (path(1:1) .eq. list_sep) then
1930 working_path_p => root_p
1931 working_path = path(2:)
1933 working_path_p => relative_p
1936 ! Loop over each field in the path
1937 do while (working_path .ne. ' ')
1938 ! Get the first list in the working path
1939 call find_head(working_path, this_list, rest)
1940 ! If the first list is empty, then the 'rest' should hold the
1941 ! final field in the path
1942 if (this_list .eq. ' ') then
1946 ! Strip off trailing list separators
1947 i = len_trim(this_list)
1948 do while (i .gt. 0 .and. this_list(i:i) .eq. list_sep)
1949 this_list(i:i) = ' '
1952 ! Find a pointer to this field in the working list
1953 this_list_p => find_field(this_list, working_path_p)
1955 if (.not. associated(this_list_p)) then
1957 ! Create the list if so requested
1958 this_list_p => make_list(working_path_p, this_list)
1959 if (.not. associated(this_list_p)) then
1964 ! Otherwise, return an error
1970 ! Make sure that the field found is a list, and if so, proceed to
1971 ! the next field in the path, otherwise, return an error
1972 if (this_list_p%field_type .eq. list_type) then
1973 working_path_p => this_list_p
1980 list_p => working_path_p
1983 end function find_list
1985 !> @brief Change the current list. Return true on success, false otherwise
1987 !> This function changes the currect list to correspond to the list named name.
1988 !! If the first character of name is the list separator (/) then the list will
1989 !! search for "name" starting from the root of the field tree. Otherwise it
1990 !! will search for name starting from the current list.
1991 !! @return A flag to indicate operation success, true = no errors
1992 function fm_change_list(name) &
1995 character(len=*), intent(in) :: name !< name of a list to change to
1997 type (field_def), pointer, save :: temp_p
1998 ! Initialize the field manager if needed
1999 if (.not. module_is_initialized) then
2000 call initialize_module_variables
2002 ! Find the list if path is not empty
2003 temp_p => find_list(name, current_list_p, .false.)
2005 if (associated(temp_p)) then
2006 current_list_p => temp_p
2012 end function fm_change_list
2014 !> @brief Change the root list
2016 !> This function changes the root of the field tree to correspond to the
2017 !! field named name. An example of a use of this would be if code is
2018 !! interested in a subset of fields with a common base. This common base
2019 !! could be set using fm_change_root and fields could be referenced using
2022 !! This function should be used in conjunction with fm_return_root.
2023 !! @return A flag to indicate operation success, true = no errors
2024 function fm_change_root(name) &
2027 character(len=*), intent(in) :: name !< name of the field which the user wishes to become the root.
2029 type (field_def), pointer, save :: temp_list_p
2031 ! Initialize the field manager if needed
2032 if (.not. module_is_initialized) then
2033 call initialize_module_variables
2036 ! Must supply a field field name
2037 if (name .eq. ' ') then
2041 ! Get a pointer to the list
2042 temp_list_p => find_list(name, current_list_p, .false.)
2044 if (associated(temp_list_p)) then
2045 ! restore the saved root values if we've already changed root
2046 if (save_root_name .ne. ' ') then
2047 root_p%name = save_root_name
2048 root_p%parent => save_root_parent_p
2050 ! set the pointer for the new root field
2051 root_p => temp_list_p
2052 ! save the new root field's name and parent
2053 save_root_name = root_p%name
2054 save_root_parent_p => root_p%parent
2055 ! set the new root name and parent fields to appropriate values
2057 nullify(root_p%parent)
2058 ! set the current list to the new root as it likely is not
2059 ! going to be meaningful anymore
2060 current_list_p => root_p
2063 ! Couldn't find the list
2067 end function fm_change_root
2069 !> @brief A function to list properties associated with a field.
2071 !> This function writes the contents of the field named "name" to stdout.
2072 !! If recursive is present and .true., then this function writes out the
2073 !! contents of any subfields associated with the field named "name".
2074 !! @return A flag to indicate operation success, true = no errors
2075 logical function fm_dump_list(name, recursive, unit) result (success)
2076 character(len=*), intent(in) :: name !< The name of the field for which output is requested.
2077 logical, intent(in), optional :: recursive !< If present and .true., then a recursive listing of
2078 !! fields will be performed.
2079 integer, intent(in), optional :: unit !< file to print to
2081 logical :: recursive_t
2082 type (field_def), pointer, save :: temp_list_p
2085 if (present(unit)) then
2091 recursive_t = .false.
2092 if (present(recursive)) recursive_t = recursive
2093 if (.not. module_is_initialized) call initialize_module_variables()
2095 if (name .eq. ' ') then
2096 ! If list is empty, then dump the current list
2097 temp_list_p => current_list_p
2100 ! Get a pointer to the list
2101 temp_list_p => find_list(name, current_list_p, .false.)
2102 if (associated(temp_list_p)) then
2110 success = dump_list(temp_list_p, recursive_t, 0, out_unit)
2112 end function fm_dump_list
2114 !> @brief A function to test whether a named field exists.
2116 !> This function determines is a field exists, relative to the current list,
2117 !! and returns true if the list exists, false otherwise.
2118 !! @return A flag to indicate operation success, true = no errors
2119 function fm_exists(name) &
2122 character(len=*), intent(in) :: name !< The name of the field that is being queried
2124 type (field_def), pointer, save :: dummy_p
2125 ! Initialize the field manager if needed
2126 if (.not. module_is_initialized) then
2127 call initialize_module_variables
2129 ! Determine whether the field exists
2130 dummy_p => get_field(name, current_list_p)
2131 success = associated(dummy_p)
2133 end function fm_exists
2135 !> @brief A function to return the index of a named field.
2137 !> Returns the index for name, returns the parameter NO_FIELD if it does not
2138 !! exist. If the first character of the named field is the list peparator,
2139 !! then the named field will be relative to the root of the field tree.
2140 !! Otherwise the named field will be relative to the current list.
2141 !> @returns index of the named field if it exists, otherwise the parameter NO_FIELD
2142 function fm_get_index(name) &
2145 character(len=*), intent(in) :: name !< The name of a field that the user wishes to get an index for
2147 type (field_def), pointer, save :: temp_field_p
2151 ! Initialize the field manager if needed
2152 if (.not. module_is_initialized) then
2153 call initialize_module_variables
2155 ! Must supply a field field name
2156 if (name .eq. ' ') then
2160 ! Get a pointer to the field
2161 temp_field_p => get_field(name, current_list_p)
2162 if (associated(temp_field_p)) then
2164 index = temp_field_p%index
2169 end function fm_get_index
2171 !> @brief A function to return the full path of the current list.
2173 !> This function returns the full path for the current list. A blank
2174 !! path indicates an error condition has occurred.
2175 !> @returns The path corresponding to the current list
2176 function fm_get_current_list() &
2178 character(len=fm_path_name_len) :: path
2180 type (field_def), pointer, save :: temp_list_p
2181 ! Initialize the field manager if needed
2182 if (.not. module_is_initialized) then
2183 call initialize_module_variables
2185 ! Set a pointer to the current list and proceed
2186 ! up the tree, filling in the name as we go
2187 temp_list_p => current_list_p
2190 do while (associated(temp_list_p))
2191 ! Check whether we are at the root field--it is the
2192 ! only field with a blank name
2193 if (temp_list_p%name .eq. ' ') then
2196 ! Append the name to the path
2197 path = list_sep // trim(temp_list_p%name) // path
2198 ! Point to the next field
2199 temp_list_p => temp_list_p%parent
2202 if (.not. associated(temp_list_p)) then
2203 ! The pointer is not associated, indicating an error has
2204 ! occurred, so set the path accordingly
2206 elseif (path .eq. ' ') then
2207 ! If path is empty, then the current list must be root,
2208 ! so set path accordingly
2212 end function fm_get_current_list
2214 !> @brief A function to return how many elements are contained within the named
2217 !> This function returns the list or entry length for the named list or entry.
2218 !! If the named field or entry does not exist, a value of 0 is returned.
2219 !> @returns The number of elements that the field name has.
2220 function fm_get_length(name) &
2223 character(len=*), intent(in) :: name !< The name of a list or entry that the user wishes to get the length of
2225 type (field_def), pointer, save :: temp_field_p
2229 ! Initialize the field manager if needed
2230 if (.not. module_is_initialized) then
2231 call initialize_module_variables
2233 ! Must supply a field name
2234 if (name .eq. ' ') then
2238 ! Get a pointer to the field
2239 temp_field_p => get_field(name, current_list_p)
2241 if (associated(temp_field_p)) then
2242 ! Set the field length
2243 if (temp_field_p%field_type .eq. list_type) then
2244 length = temp_field_p%length
2246 length = temp_field_p%max_index
2252 end function fm_get_length
2254 !> @brief A function to return the type of the named field
2256 !> This function returns the type of the field for name.
2257 !! This indicates whether the named field is a "list" (has children fields),
2258 !! or has values of type "integer", "real", "logical" or "string".
2259 !! If it does not exist it returns a blank string.
2260 !> @returns A string containing the type of the named field
2261 function fm_get_type(name) &
2262 result (name_field_type)
2263 character(len=8) :: name_field_type
2264 character(len=*), intent(in) :: name !< The name of a field that the user wishes to find the type of
2266 type (field_def), pointer, save :: temp_field_p
2270 ! Initialize the field manager if needed
2271 if (.not. module_is_initialized) then
2272 call initialize_module_variables
2274 ! Must supply a field name
2275 if (name .eq. ' ') then
2276 name_field_type = ' '
2279 ! Get a pointer to the field
2280 temp_field_p => get_field(name, current_list_p)
2282 if (associated(temp_field_p)) then
2283 ! Set the field type
2284 name_field_type = field_type_name(temp_field_p%field_type)
2286 name_field_type = ' '
2289 end function fm_get_type
2291 !> @returns A flag to indicate whether the function operated with (false) or without
2293 function fm_get_value_integer(name, get_ival, index) &
2296 character(len=*), intent(in) :: name !< The name of a field that the user wishes to get a value for.
2297 integer, intent(out) :: get_ival !< The value associated with the named field.
2298 integer, intent(in), optional :: index !< An optional index to retrieve a single value from an array.
2301 type (field_def), pointer, save :: temp_field_p
2305 ! Initialize the field manager if needed
2306 if (.not. module_is_initialized) then
2307 call initialize_module_variables
2309 ! Must supply a field field name
2310 if (name .eq. ' ') then
2315 ! Set index to retrieve
2316 if (present(index)) then
2321 ! Get a pointer to the field
2322 temp_field_p => get_field(name, current_list_p)
2324 if (associated(temp_field_p)) then
2325 ! check that the field is the correct type
2326 if (temp_field_p%field_type .eq. integer_type) then
2327 if (index_t .lt. 1 .or. index_t .gt. temp_field_p%max_index) then
2328 ! Index is not positive or index too large
2333 get_ival = temp_field_p%i_value(index_t)
2337 ! Field not corrcet type
2346 end function fm_get_value_integer
2348 !> @returns A flag to indicate whether the function operated with (false) or without
2350 function fm_get_value_logical(name, get_lval, index) &
2353 character(len=*), intent(in) :: name !< The name of a field that the user wishes to get a value for.
2354 logical, intent(out) :: get_lval !< The value associated with the named field
2355 integer, intent(in), optional :: index !< An optional index to retrieve a single value from an array.
2358 type (field_def), pointer, save :: temp_field_p
2362 ! Initialize the field manager if needed
2363 if (.not. module_is_initialized) then
2364 call initialize_module_variables
2366 ! Must supply a field field name
2367 if (name .eq. ' ') then
2372 ! Set index to retrieve
2373 if (present(index)) then
2378 ! Get a pointer to the field
2379 temp_field_p => get_field(name, current_list_p)
2381 if (associated(temp_field_p)) then
2382 ! check that the field is the correct type
2383 if (temp_field_p%field_type .eq. logical_type) then
2385 if (index_t .lt. 1 .or. index_t .gt. temp_field_p%max_index) then
2386 ! Index is not positive or too large
2391 get_lval = temp_field_p%l_value(index_t)
2395 ! Field not correct type
2404 end function fm_get_value_logical
2406 !> @returns A flag to indicate whether the function operated with (false) or without
2408 function fm_get_value_string(name, get_sval, index) &
2411 character(len=*), intent(in) :: name !< The name of a field that the user wishes to get a value for.
2412 character(len=*), intent(out) :: get_sval !< The value associated with the named field
2413 integer, intent(in), optional :: index !< An optional index to retrieve a single value from an array.
2416 type (field_def), pointer, save :: temp_field_p
2420 ! Initialize the field manager if needed
2421 if (.not. module_is_initialized) then
2422 call initialize_module_variables
2424 ! Must supply a field field name
2425 if (name .eq. ' ') then
2430 ! Set index to retrieve
2431 if (present(index)) then
2436 ! Get a pointer to the field
2437 temp_field_p => get_field(name, current_list_p)
2439 if (associated(temp_field_p)) then
2440 ! check that the field is the correct type
2441 if (temp_field_p%field_type .eq. string_type) then
2442 if (index_t .lt. 1 .or. index_t .gt. temp_field_p%max_index) then
2443 ! Index is not positive or is too large
2448 get_sval = temp_field_p%s_value(index_t)
2452 ! Field not correct type
2461 end function fm_get_value_string
2463 !> Iterates through the given list
2464 !> @returns A flag to indicate whether the function operated with (FALSE)
2465 !! or without (TRUE) errors
2466 function fm_loop_over_list_old(list, name, field_type, index) &
2469 character(len=*), intent(in) :: list !< Name of a list to loop over
2470 character(len=*), intent(out) :: name !< name of a field from list
2471 character(len=fm_type_name_len), intent(out) :: field_type !< type of a list entry
2472 integer, intent(out) :: index !< index of the field within the list
2477 ! Initialize the field manager if needed
2478 if (.not. module_is_initialized) then
2479 call initialize_module_variables
2482 if (list .eq. loop_list .and. associated(loop_list_p)) then
2483 ! We've already started this loop, so continue on
2484 loop_list_p => loop_list_p%next
2485 success = set_list_stuff()
2486 elseif (list .eq. ' ') then
2487 ! If list is empty, then loop over the current list
2489 loop_list_p => current_list_p%first_field
2490 success = set_list_stuff()
2492 ! Get a pointer to the list
2494 loop_list_p => find_list(loop_list, current_list_p, .false.)
2495 if (associated(loop_list_p)) then
2496 loop_list_p => loop_list_p%first_field
2497 success = set_list_stuff()
2507 !> If the the pointer matches to the right list,
2508 !! extract the field information. Used in fm_loop_over_list
2509 !> @returns A flag to indicate whether the function operated with (FALSE)
2510 !! or without (TRUE) errors
2511 function set_list_stuff() &
2516 if (associated(loop_list_p)) then
2517 name = loop_list_p%name
2518 field_type = field_type_name(loop_list_p%field_type)
2519 index = loop_list_p%index
2529 end function set_list_stuff
2531 end function fm_loop_over_list_old
2533 !> given a name of the list, prepares an iterator over the list content.
2534 !! If the name of the given list is blank, then the current list is used
2535 subroutine fm_init_loop(loop_list, iter)
2536 character(len=*) , intent(in) :: loop_list !< name of the list to iterate over
2537 type(fm_list_iter_type), intent(out) :: iter !< loop iterator
2539 if (.not.module_is_initialized) call initialize_module_variables
2541 if (loop_list==' ') then ! looping over current list
2542 iter%ptr => current_list_p%first_field
2544 iter%ptr => find_list(loop_list,current_list_p,.false.)
2545 if (associated(iter%ptr)) iter%ptr => iter%ptr%first_field
2547 end subroutine fm_init_loop
2549 !> given a list iterator, returns information about curren list element
2550 !! and advances the iterator to the next list element. At the end of the
2551 !! list, returns FALSE
2552 function fm_loop_over_list_new(iter, name, field_type, index) &
2553 result (success) ; logical success
2554 type (fm_list_iter_type), intent(inout) :: iter !< list iterator
2555 character(len=*), intent(out) :: name !< name of the current list item
2556 character(len=*), intent(out) :: field_type !< type of the field
2557 integer , intent(out) :: index !< index in the list
2559 if (.not.module_is_initialized) call initialize_module_variables
2560 if (associated(iter%ptr)) then
2561 name = iter%ptr%name
2562 field_type = field_type_name(iter%ptr%field_type)
2563 index = iter%ptr%index
2565 iter%ptr => iter%ptr%next
2572 end function fm_loop_over_list_new
2574 !> @brief A function to create a new list
2576 !> Allocate and initialize a new list and return the index of the list. If an
2577 !! error occurs return the parameter NO_FIELD.
2578 !> @return integer index of the newly created list
2579 function fm_new_list(name, create, keep) &
2582 character(len=*), intent(in) :: name !< Name of a list that user wishes to create
2583 logical, intent(in), optional :: create !< If present and true, create the list if it does not exist
2584 logical, intent(in), optional :: keep !< If present and true, make this list the current list
2588 character(len=fm_path_name_len) :: path
2589 character(len=fm_field_name_len) :: base
2590 type (field_def), pointer, save :: temp_list_p
2594 ! Initialize the field manager if needed
2595 if (.not. module_is_initialized) then
2596 call initialize_module_variables
2598 ! Must supply a field list name
2599 if (name .eq. ' ') then
2603 ! Check for optional arguments
2604 if (present(create)) then
2610 if (present(keep)) then
2615 ! Get a pointer to the parent list
2616 call find_base(name, path, base)
2618 temp_list_p => find_list(path, current_list_p, create_t)
2620 if (associated(temp_list_p)) then
2622 temp_list_p => make_list(temp_list_p, base)
2623 if (associated(temp_list_p)) then
2624 ! Make this list the current list, if requested
2626 current_list_p => temp_list_p
2628 index = temp_list_p%index
2636 end function fm_new_list
2638 !> @brief Assigns a given value to a given field
2639 !> @returns An index for the named field
2640 function fm_new_value_integer(name, new_ival, create, index, append) &
2641 result (field_index)
2642 integer :: field_index
2643 character(len=*), intent(in) :: name !< The name of a field that the user wishes to create
2645 integer, intent(in) :: new_ival !< The value that the user wishes to apply to the
2647 logical, intent(in), optional :: create !< If present and .true., then a value for this
2648 !! field will be created.
2649 integer, intent(in), optional :: index !< The index to an array of values that the user
2650 !! wishes to apply a new value.
2651 logical, intent(in), optional :: append !< If present and .true., then append the value to
2652 !! an array of the present values. If present and .true., then index cannot be greater than 0.
2657 integer, pointer, dimension(:) :: temp_i_value
2658 character(len=fm_path_name_len) :: path
2659 character(len=fm_field_name_len) :: base
2660 type (field_def), pointer, save :: temp_list_p
2661 type (field_def), pointer, save :: temp_field_p
2665 ! Initialize the field manager if needed
2666 if (.not. module_is_initialized) then
2667 call initialize_module_variables
2669 ! Must supply a field name
2670 if (name .eq. ' ') then
2671 field_index = NO_FIELD
2674 ! Check for optional arguments
2675 if (present(create)) then
2680 ! Check that append is not true and index non-positive
2681 if (present(index) .and. present(append)) then
2682 if (append .and. index .gt. 0) then
2683 field_index = NO_FIELD
2687 ! Set index to define
2688 if (present(index)) then
2690 if (index_t .lt. 0) then
2692 field_index = NO_FIELD
2698 ! Get a pointer to the parent list
2699 call find_base(name, path, base)
2700 temp_list_p => find_list(path, current_list_p, create_t)
2702 if (associated(temp_list_p)) then
2703 temp_field_p => find_field(base, temp_list_p)
2704 if (.not. associated(temp_field_p)) then
2705 ! Create the field if it doesn't exist
2706 temp_field_p => create_field(temp_list_p, base)
2708 if (associated(temp_field_p)) then
2709 ! Check if the field_type is the same as previously
2710 ! If not then reset max_index to 0
2711 if (temp_field_p%field_type == real_type ) then
2712 ! promote integer input to real
2713 ! all real field values are stored as r8_kind
2714 field_index = fm_new_value(name, real(new_ival,r8_kind), create, index, append)
2716 else if (temp_field_p%field_type /= integer_type ) then
2717 ! slm: why would we reset index? Is it not an error to have a "list" defined
2718 ! with different types in more than one place?
2719 temp_field_p%max_index = 0
2722 temp_field_p%field_type = integer_type
2723 ! Set the index if appending
2724 if (present(append)) then
2726 index_t = temp_field_p%max_index + 1
2730 if (index_t .gt. temp_field_p%max_index + 1) then
2732 field_index = NO_FIELD
2735 elseif (index_t .eq. 0 .and. &
2736 temp_field_p%max_index .gt. 0) then
2737 ! Can't set non-null field to null
2738 field_index = NO_FIELD
2741 elseif (.not. allocated(temp_field_p%i_value) .and. &
2742 index_t .gt. 0) then
2743 ! Array undefined, so allocate the array
2744 allocate(temp_field_p%i_value(1))
2745 temp_field_p%max_index = 1
2746 temp_field_p%array_dim = 1
2747 elseif (index_t .gt. temp_field_p%array_dim) then
2748 ! Array is too small, so allocate new array and copy over
2750 temp_field_p%array_dim = temp_field_p%array_dim + array_increment
2751 allocate (temp_i_value(temp_field_p%array_dim))
2752 do i = 1, temp_field_p%max_index
2753 temp_i_value(i) = temp_field_p%i_value(i)
2755 if (allocated(temp_field_p%i_value)) deallocate(temp_field_p%i_value)
2756 temp_field_p%i_value = temp_i_value
2757 temp_field_p%max_index = index_t
2759 ! Assign the value and set the field_index for return
2760 ! for non-null fields (index_t > 0)
2761 if (index_t .gt. 0) then
2762 temp_field_p%i_value(index_t) = new_ival
2763 if (index_t .gt. temp_field_p%max_index) then
2764 temp_field_p%max_index = index_t
2767 field_index = temp_field_p%index
2770 field_index = NO_FIELD
2773 field_index = NO_FIELD
2776 end function fm_new_value_integer
2778 !> @brief Assigns a given value to a given field
2779 !> @returns An index for the named field
2780 function fm_new_value_logical(name, new_lval, create, index, append) &
2781 result (field_index)
2782 integer :: field_index
2783 character(len=*), intent(in) :: name !< The name of a field that the user wishes to create
2785 logical, intent(in) :: new_lval !< The value that the user wishes to apply to the
2787 logical, intent(in), optional :: create !< If present and .true., then a value for this
2788 !! field will be created.
2789 integer, intent(in), optional :: index !< The index to an array of values that the user
2790 !! wishes to apply a new value.
2791 logical, intent(in), optional :: append !< If present and .true., then append the value to
2792 !! an array of the present values. If present and .true., then index cannot be greater than 0.
2794 character(len=fm_path_name_len) :: path
2795 character(len=fm_field_name_len) :: base
2799 logical, dimension(:), pointer :: temp_l_value
2800 type (field_def), pointer, save :: temp_list_p
2801 type (field_def), pointer, save :: temp_field_p
2805 ! Initialize the field manager if needed
2806 if (.not. module_is_initialized) then
2807 call initialize_module_variables
2809 ! Must supply a field name
2810 if (name .eq. ' ') then
2811 field_index = NO_FIELD
2814 ! Check for optional arguments
2815 if (present(create)) then
2820 ! Check that append is not true and index greater than 0
2821 if (present(index) .and. present(append)) then
2822 if (append .and. index .gt. 0) then
2823 field_index = NO_FIELD
2827 ! Set index to define
2828 if (present(index)) then
2830 if (index_t .lt. 0) then
2832 field_index = NO_FIELD
2838 ! Get a pointer to the parent list
2839 call find_base(name, path, base)
2840 temp_list_p => find_list(path, current_list_p, create_t)
2842 if (associated(temp_list_p)) then
2843 temp_field_p => find_field(base, temp_list_p)
2844 if (.not. associated(temp_field_p)) then
2845 ! Create the field if it doesn't exist
2846 temp_field_p => create_field(temp_list_p, base)
2848 if (associated(temp_field_p)) then
2849 ! Check if the field_type is the same as previously
2850 ! If not then reset max_index to 0
2851 if (temp_field_p%field_type /= logical_type ) then
2852 temp_field_p%max_index = 0
2855 temp_field_p%field_type = logical_type
2856 ! Set the index if appending
2857 if (present(append)) then
2859 index_t = temp_field_p%max_index + 1
2863 if (index_t .gt. temp_field_p%max_index + 1) then
2865 field_index = NO_FIELD
2868 elseif (index_t .eq. 0 .and. &
2869 temp_field_p%max_index .gt. 0) then
2870 ! Can't set non-null field to null
2871 field_index = NO_FIELD
2874 elseif (.not. allocated(temp_field_p%l_value) .and. &
2875 index_t .gt. 0) then
2876 ! Array undefined, so allocate the array
2877 allocate(temp_field_p%l_value(1))
2878 temp_field_p%max_index = 1
2879 temp_field_p%array_dim = 1
2881 elseif (index_t .gt. temp_field_p%array_dim) then
2882 ! Array is too small, so allocate new array and copy over
2884 temp_field_p%array_dim = temp_field_p%array_dim + array_increment
2885 allocate (temp_l_value(temp_field_p%array_dim))
2886 do i = 1, temp_field_p%max_index
2887 temp_l_value(i) = temp_field_p%l_value(i)
2889 if (allocated(temp_field_p%l_value)) deallocate(temp_field_p%l_value)
2890 temp_field_p%l_value = temp_l_value
2891 temp_field_p%max_index = index_t
2894 ! Assign the value and set the field_index for return
2895 ! for non-null fields (index_t > 0)
2896 if (index_t .gt. 0) then
2897 temp_field_p%l_value(index_t) = new_lval
2898 if (index_t .gt. temp_field_p%max_index) then
2899 temp_field_p%max_index = index_t
2902 field_index = temp_field_p%index
2904 field_index = NO_FIELD
2907 field_index = NO_FIELD
2910 end function fm_new_value_logical
2912 !> @brief Assigns a given value to a given field
2913 !> @returns An index for the named field
2914 function fm_new_value_string(name, new_sval, create, index, append) &
2915 result (field_index)
2916 integer :: field_index
2917 character(len=*), intent(in) :: name !< The name of a field that the user wishes to create
2919 character(len=*), intent(in) :: new_sval !< The value that the user wishes to apply to the
2921 logical, intent(in), optional :: create !< If present and .true., then a value for this
2922 !! field will be created.
2923 integer, intent(in), optional :: index !< The index to an array of values that the user
2924 !! wishes to apply a new value.
2925 logical, intent(in), optional :: append !< If present and .true., then append the value to
2927 character(len=fm_string_len), dimension(:), pointer :: temp_s_value
2928 character(len=fm_path_name_len) :: path
2929 character(len=fm_field_name_len) :: base
2933 type (field_def), save, pointer :: temp_list_p
2934 type (field_def), save, pointer :: temp_field_p
2938 ! Initialize the field manager if needed
2939 if (.not. module_is_initialized) then
2940 call initialize_module_variables
2942 ! Must supply a field name
2943 if (name .eq. ' ') then
2944 field_index = NO_FIELD
2947 ! Check for optional arguments
2948 if (present(create)) then
2953 ! Check that append is not true and index greater than 0
2954 if (present(index) .and. present(append)) then
2955 if (append .and. index .gt. 0) then
2956 field_index = NO_FIELD
2960 ! Set index to define
2961 if (present(index)) then
2963 if (index_t .lt. 0) then
2965 field_index = NO_FIELD
2971 ! Get a pointer to the parent list
2972 call find_base(name, path, base)
2973 temp_list_p => find_list(path, current_list_p, create_t)
2975 if (associated(temp_list_p)) then
2976 temp_field_p => find_field(base, temp_list_p)
2977 if (.not. associated(temp_field_p)) then
2978 ! Create the field if it doesn't exist
2979 temp_field_p => create_field(temp_list_p, base)
2981 if (associated(temp_field_p)) then
2982 ! Check if the field_type is the same as previously
2983 ! If not then reset max_index to 0
2984 if (temp_field_p%field_type /= string_type ) then
2985 temp_field_p%max_index = 0
2988 temp_field_p%field_type = string_type
2989 ! Set the index if appending
2990 if (present(append)) then
2992 index_t = temp_field_p%max_index + 1
2996 if (index_t .gt. temp_field_p%max_index + 1) then
2998 field_index = NO_FIELD
3001 elseif (index_t .eq. 0 .and. &
3002 temp_field_p%max_index .gt. 0) then
3003 ! Can't set non-null field to null
3004 field_index = NO_FIELD
3007 elseif (.not.allocated(temp_field_p%s_value) .and. &
3008 index_t .gt. 0) then
3009 ! Array undefined, so allocate the array
3010 allocate(temp_field_p%s_value(1))
3011 temp_field_p%max_index = 1
3012 temp_field_p%array_dim = 1
3014 elseif (index_t .gt. temp_field_p%array_dim) then
3015 ! Array is too small, so allocate new array and copy over
3017 temp_field_p%array_dim = temp_field_p%array_dim + array_increment
3018 allocate (temp_s_value(temp_field_p%array_dim))
3019 do i = 1, temp_field_p%max_index
3020 temp_s_value(i) = temp_field_p%s_value(i)
3022 if (allocated(temp_field_p%s_value)) deallocate(temp_field_p%s_value)
3023 temp_field_p%s_value = temp_s_value
3024 temp_field_p%max_index = index_t
3027 ! Assign the value and set the field_index for return
3028 ! for non-null fields (index_t > 0)
3029 if (index_t .gt. 0) then
3030 temp_field_p%s_value(index_t) = new_sval
3031 if (index_t .gt. temp_field_p%max_index) then
3032 temp_field_p%max_index = index_t
3035 field_index = temp_field_p%index
3037 ! Error in making the field
3038 field_index = NO_FIELD
3041 ! Error following the path
3042 field_index = NO_FIELD
3045 end function fm_new_value_string
3048 !> Resets the loop variable. For use in conjunction with fm_loop_over_list.
3049 subroutine fm_reset_loop
3050 ! Initialize the field manager if needed
3051 if (.not. module_is_initialized) then
3052 call initialize_module_variables
3054 ! Reset the variables
3056 nullify(loop_list_p)
3058 end subroutine fm_reset_loop
3060 !> Return the root list to the value at initialization.
3062 !> For use in conjunction with fm_change_root.
3064 !! Users should use this routine before leaving their routine if they
3065 !! previously used fm_change_root.
3066 subroutine fm_return_root
3067 ! Initialize the field manager if needed
3068 if (.not. module_is_initialized) then
3069 call initialize_module_variables
3071 ! restore the saved values to the current root
3072 root_p%name = save_root_name
3073 root_p%parent => save_root_parent_p
3074 ! set the pointer to the original root field
3076 ! reset the save root name and parent variables
3077 save_root_name = ' '
3078 nullify(save_root_parent_p)
3080 end subroutine fm_return_root
3082 !> Return a pointer to the field if it exists relative to this_list_p,
3084 !! @returns A pointer to the field name
3085 function get_field(name, this_list_p) &
3087 type (field_def), pointer :: list_p
3088 character(len=*), intent(in) :: name !< The name of a list that the user wishes to get information for
3089 type (field_def), pointer :: this_list_p !< A pointer to a list that serves as the base point
3090 !! for searching for name
3092 character(len=fm_path_name_len) :: path
3093 character(len=fm_field_name_len) :: base
3094 type (field_def), pointer, save :: temp_p
3097 ! Get the path and base for name
3098 call find_base(name, path, base)
3099 ! Find the list if path is not empty
3100 if (path .ne. ' ') then
3101 temp_p => find_list(path, this_list_p, .false.)
3102 if (associated(temp_p)) then
3103 list_p => find_field(base, temp_p)
3108 list_p => find_field(base, this_list_p)
3111 end function get_field
3113 !> This function allows a user to rename a field without modifying the
3114 !! contents of the field.
3116 !> Function to modify the name of a field.
3117 !! Should be used with caution.
3118 !! @returns A flag to indicate whether the function operated with (FALSE) or
3119 !! without (TRUE) errors.
3120 function fm_modify_name(oldname, newname) &
3123 character(len=*), intent(in) :: oldname !< The name of a field that the user wishes to change
3125 character(len=*), intent(in) :: newname !< The name that the user wishes to change the name of
3128 character(len=fm_path_name_len) :: path
3129 character(len=fm_field_name_len) :: base
3130 type (field_def), pointer, save :: list_p
3131 type (field_def), pointer, save :: temp_p
3132 ! Get the path and base for name
3133 call find_base(oldname, path, base)
3134 ! Find the list if path is not empty
3136 if (path .ne. ' ') then
3137 temp_p => find_list(path, current_list_p, .false.)
3138 if (associated(temp_p)) then
3139 list_p => find_field(base, temp_p)
3140 if (associated(list_p)) then
3141 list_p%name = newname
3148 list_p => find_field(base, current_list_p)
3149 if (associated(list_p)) then
3150 list_p%name = newname
3155 end function fm_modify_name
3157 !> A function to initialize the values of the pointers. This will remove
3158 !! all fields and reset the field tree to only the root field.
3159 subroutine initialize_module_variables
3160 ! Initialize the root field
3161 integer :: io, ierr !< Error codes when reading the namelist
3162 integer :: logunit !< Unit number for the log file
3164 if (.not. module_is_initialized) then
3166 read (input_nml_file, nml=field_manager_nml, iostat=io)
3167 ierr = check_nml_error(io,"field_manager_nml")
3170 if (mpp_pe() == mpp_root_pe()) write (logunit, nml=field_manager_nml)
3174 field_type_name(integer_type) = 'integer'
3175 field_type_name(list_type) = 'list'
3176 field_type_name(logical_type) = 'logical'
3177 field_type_name(real_type) = 'real'
3178 field_type_name(string_type) = 'string'
3182 root%parent => root_p
3184 root%field_type = list_type
3187 nullify(root%first_field)
3188 nullify(root%last_field)
3191 if (allocated(root%i_value)) deallocate(root%i_value)
3192 if (allocated(root%l_value)) deallocate(root%l_value)
3193 if (allocated(root%r_value)) deallocate(root%r_value)
3194 if (allocated(root%s_value)) deallocate(root%s_value)
3199 current_list_p => root
3201 nullify(loop_list_p)
3204 nullify(save_root_parent_p)
3205 save_root_name = ' '
3207 module_is_initialized = .true.
3211 end subroutine initialize_module_variables
3213 !> This function creates a new field and returns a pointer to that field.
3215 !> Allocate and initialize a new list in this_list_p list.
3216 !! @return a pointer to the list on success, or a null pointer
3218 function make_list(this_list_p, name) &
3220 type (field_def), pointer :: list_p
3221 type (field_def), pointer :: this_list_p !< Base of a list that the user wishes to add a list to
3222 character(len=*), intent(in) :: name !< name of a list that the user wishes to create
3224 type (field_def), pointer, save :: dummy_p
3228 ! Check to see whether there is already a list with
3229 ! this name, and if so, return an error as list names
3231 dummy_p => find_field(name, this_list_p )
3232 if (associated(dummy_p)) then
3233 ! This list is already specified, return an error
3237 ! Create a field for the new list
3239 list_p => create_field(this_list_p, name)
3240 if (.not. associated(list_p)) then
3244 ! Initialize the new list
3246 list_p%field_type = list_type
3247 if (allocated(list_p%i_value)) deallocate(list_p%i_value)
3248 if (allocated(list_p%l_value)) deallocate(list_p%l_value)
3249 if (allocated(list_p%r_value)) deallocate(list_p%r_value)
3250 if (allocated(list_p%s_value)) deallocate(list_p%s_value)
3252 end function make_list
3254 !> This is a function that provides the capability to return parameters
3255 !! associated with a field in a pair of strings.
3257 !> Given a name return a list of method names and control strings.
3258 !! This function should return strings similar to those in the field
3259 !! table if a comma delimited format is being used.
3260 !> @return A flag to indicate whether the function operated with (FALSE) or
3261 !! without (TRUE) errors
3262 function fm_query_method(name, method_name, method_control) &
3265 character(len=*), intent(in) :: name !< name of a list that the user wishes to change to
3266 character(len=*), intent(out) :: method_name !< name of a parameter associated with the named field
3267 character(len=*), intent(out) :: method_control !< value of parameters associated with the named field
3269 character(len=fm_path_name_len) :: path
3270 character(len=fm_path_name_len) :: base
3271 character(len=fm_path_name_len) :: name_loc
3272 logical :: recursive_t
3273 type (field_def), pointer, save :: temp_list_p
3274 type (field_def), pointer, save :: temp_value_p
3275 type (field_def), pointer, save :: this_field_p
3280 recursive_t = .true.
3282 method_control = " "
3283 ! Initialize the field manager if needed
3284 if (.not. module_is_initialized) call initialize_module_variables
3285 name_loc = lowercase(name)
3286 call find_base(name_loc, path, base)
3288 temp_list_p => find_list(name_loc, current_list_p, .false.)
3290 if (associated(temp_list_p)) then
3291 ! Find the entry values for the list.
3292 success = query_method(temp_list_p, recursive_t, base, method_name, method_control)
3294 ! This is not a list but it may be a parameter with a value
3295 ! If so put the parameter value in method_name.
3297 temp_value_p => find_list(path, current_list_p, .false.)
3298 if (associated(temp_value_p)) then
3299 ! Find the entry values for this item.
3300 this_field_p => temp_value_p%first_field
3302 do while (associated(this_field_p))
3303 if ( this_field_p%name == base ) then
3304 method_name = this_field_p%s_value(1)
3311 this_field_p => this_field_p%next
3315 ! Error following the path
3320 end function fm_query_method
3322 !> A private function that can recursively recover values for parameters
3323 !! associated with a field.
3324 !> @return A flag to indicate whether the function operated with (FALSE) or
3325 !! without (TRUE) errors
3326 recursive function query_method(list_p, recursive, name, method_name, method_control) &
3329 type (field_def), pointer :: list_p !< A pointer to the field that is of interest
3330 logical, intent(in) :: recursive !< A flag to enable recursive searching if true
3331 character(len=*), intent(in) :: name !< name of a list that the user wishes to change to
3332 character(len=*), intent(out) :: method_name !< name of a parameter associated with the named field
3333 character(len=*), intent(out) :: method_control !< value of parameters associated with the named field
3336 character(len=64) :: scratch
3337 type (field_def), pointer :: this_field_p
3342 ! Check for a valid list
3343 if (.not. associated(list_p) .or. list_p%field_type .ne. list_type) then
3347 ! set the default return value
3350 this_field_p => list_p%first_field
3352 do while (associated(this_field_p))
3353 select case(this_field_p%field_type)
3355 ! If this is a list, then this is the method name
3357 if (.not. query_method(this_field_p, .true., this_field_p%name, method_name, method_control)) then
3361 method_name = trim(method_name)//trim(this_field_p%name)
3366 write (scratch,*) this_field_p%i_value
3367 call concat_strings(method_control, comma//trim(this_field_p%name)//' = '//trim(adjustl(scratch)))
3370 write (scratch,'(l1)')this_field_p%l_value
3371 call concat_strings(method_control, comma//trim(this_field_p%name)//' = '//trim(adjustl(scratch)))
3374 write (scratch,*) this_field_p%r_value
3375 call concat_strings(method_control, comma//trim(this_field_p%name)//' = '//trim(adjustl(scratch)))
3378 call concat_strings(method_control, comma//trim(this_field_p%name)//' = '//trim(this_field_p%s_value(1)))
3379 do i = 2, this_field_p%max_index
3380 call concat_strings(method_control, comma//trim(this_field_p%s_value(i)))
3388 this_field_p => this_field_p%next
3392 end function query_method
3394 !> private function: appends str2 to the end of str1, with length check
3395 subroutine concat_strings(str1,str2)
3396 character(*), intent(inout) :: str1
3397 character(*), intent(in) :: str2
3399 character(64) :: n1,n2 ! for error reporting
3401 if (len_trim(str1)+len_trim(str2)>len(str1)) then
3402 write(n1,*)len(str1)
3403 write(n2,*)len_trim(str1)+len_trim(str2)
3404 call mpp_error(FATAL,'length of output string ('//trim(adjustl(n1))&
3405 //') is not enough for the result of concatenation (len='&
3406 //trim(adjustl(n2))//')')
3408 str1 = trim(str1)//trim(str2)
3409 end subroutine concat_strings
3411 !> A function that allows the user to copy a field and add a suffix to
3412 !! the name of the new field.
3414 !> Given the name of a pre-existing field and a suffix, this function
3415 !! will create a new field. The name of the new field will be that of
3416 !! the old field with a suffix supplied by the user.
3417 !! @return index of the field that has been created by the copy
3418 function fm_copy_list(list_name, suffix, create ) &
3421 character(len=*), intent(in) :: list_name !< name of a field that the user wishes to copy
3422 character(len=*), intent(in) :: suffix !< suffix that will be added to list_name when
3424 logical, intent(in), optional :: create !< flag to create new list if applicable
3426 character(len=fm_string_len), dimension(:), allocatable :: control
3427 character(len=fm_string_len), dimension(:), allocatable :: method
3428 character(len=fm_string_len) :: head
3429 character(len=fm_string_len) :: list_name_new
3430 character(len=fm_string_len) :: tail
3431 character(len=fm_string_len) :: val_str
3435 logical :: found_methods
3436 logical :: got_value
3437 logical :: recursive_t
3439 logical :: val_logical
3440 real(r8_kind) :: val_real
3441 type (field_def), pointer, save :: temp_field_p
3442 type (field_def), pointer, save :: temp_list_p
3449 list_name_new = trim(list_name)//trim(suffix)
3450 recursive_t = .true.
3451 ! Initialize the field manager if needed
3452 if (.not. module_is_initialized) then
3453 call initialize_module_variables
3456 if (list_name .eq. ' ') then
3457 ! If list is empty, then dump the current list
3458 temp_list_p => current_list_p
3461 ! Get a pointer to the list
3462 temp_list_p => find_list(list_name, current_list_p, .false.)
3463 if (associated(temp_list_p)) then
3466 ! Error following the path
3472 found_methods = fm_find_methods(trim(list_name), method, control)
3473 do n = 1, size(method)
3474 if (LEN_TRIM(method(n)) > 0 ) then
3475 index = fm_new_list(trim(list_name_new)//list_sep//method(n), create = create)
3476 call find_base(method(n), head, tail)
3477 temp_field_p => find_list(trim(list_name)//list_sep//head,temp_list_p, .false.)
3478 temp_field_p => find_field(tail,temp_field_p)
3479 select case (temp_field_p%field_type)
3481 got_value = fm_get_value( trim(list_name)//list_sep//method(n), val_int)
3482 if ( fm_new_value( trim(list_name_new)//list_sep//method(n), val_int, &
3483 create = create, append = .true.) < 0 ) &
3484 call mpp_error(FATAL, trim(error_header)//'Could not set the '//trim(method(n))//&
3485 ' for '//trim(list_name)//trim(suffix))
3488 got_value = fm_get_value( trim(list_name)//list_sep//method(n), val_logical)
3489 if ( fm_new_value( trim(list_name_new)//list_sep//method(n), val_logical, &
3490 create = create, append = .true.) < 0 ) &
3491 call mpp_error(FATAL, trim(error_header)//'Could not set the '//trim(method(n))//&
3492 ' for '//trim(list_name)//trim(suffix))
3495 got_value = fm_get_value( trim(list_name)//list_sep//method(n), val_real)
3496 if ( fm_new_value( trim(list_name_new)//list_sep//method(n), val_real, &
3497 create = create, append = .true.) < 0 ) &
3498 call mpp_error(FATAL, trim(error_header)//'Could not set the '//trim(method(n))//&
3499 ' for '//trim(list_name)//trim(suffix))
3502 got_value = fm_get_value( trim(list_name)//list_sep//method(n), val_str)
3503 if ( fm_new_value( trim(list_name_new)//list_sep//method(n), val_str, &
3504 create = create, append = .true.) < 0 ) &
3505 call mpp_error(FATAL, trim(error_header)//'Could not set the '//trim(method(n))//&
3506 ' for '//trim(list_name)//trim(suffix))
3514 end function fm_copy_list
3516 !> This function retrieves all the methods associated with a field.
3518 !> This is different from fm_query_method in that this function gets all
3519 !! the methods associated as opposed to 1 method.
3520 !! @return A flag to indicate whether the function operated with (FALSE) or
3521 !! without (TRUE) errors.
3522 function fm_find_methods(list_name, methods, control ) &
3525 character(len=*), intent(in) :: list_name !< The name of a list that the user wishes to find methods for
3526 character(len=*), intent(out), dimension(:) :: methods !< An array of the methods associated with list_name
3527 character(len=*), intent(out), dimension(:) :: control !< An array of the parameters associated with methods
3530 logical :: recursive_t
3531 type (field_def), pointer, save :: temp_list_p
3536 ! Check whether to do things recursively
3537 recursive_t = .true.
3539 if (.not. module_is_initialized) then
3540 call initialize_module_variables
3543 if (list_name .eq. ' ') then
3544 ! If list is empty, then dump the current list
3545 temp_list_p => current_list_p
3548 ! Get a pointer to the list
3549 temp_list_p => find_list(list_name, current_list_p, .false.)
3550 if (associated(temp_list_p)) then
3553 ! Error following the path
3559 success = find_method(temp_list_p, recursive_t, num_meth, methods, control)
3562 end function fm_find_methods
3564 !> Given a field list pointer this function retrieves methods and
3565 !! associated parameters for the field list.
3566 !! @returns A flag to indicate whether the function operated with (FALSE) or
3567 !! without (TRUE) errors.
3568 recursive function find_method(list_p, recursive, num_meth, method, control) &
3571 type (field_def), pointer :: list_p !< A pointer to the field of interest
3572 logical, intent(in) :: recursive !< If true, search recursively for fields
3573 integer, intent(inout) :: num_meth !< The number of methods found
3574 character(len=*), intent(out), dimension(:) :: method !< The methods associated with the field pointed to by list_p
3575 character(len=*), intent(out), dimension(:) :: control !< The control parameters for the methods found
3577 character(len=fm_path_name_len) :: scratch
3580 type (field_def), pointer, save :: this_field_p
3584 ! Check for a valid list
3585 if (.not. associated(list_p) .or. list_p%field_type .ne. list_type) then
3588 ! set the default return value
3591 this_field_p => list_p%first_field
3593 do while (associated(this_field_p))
3594 select case(this_field_p%field_type)
3596 ! If this is a list, then this is the method name
3597 if ( this_field_p%length > 1) then
3598 do n = num_meth+1, num_meth + this_field_p%length - 1
3599 write (method(n),'(a,a,a,$)') trim(method(num_meth)), &
3600 trim(this_field_p%name), list_sep
3602 write (method(num_meth),'(a,a,a,$)') trim(method(num_meth)), &
3603 trim(this_field_p%name), list_sep
3605 write (method(num_meth),'(a,a,a,$)') trim(method(num_meth)), &
3606 trim(this_field_p%name), list_sep
3608 success = find_method(this_field_p, .true., num_meth, method, control)
3611 write (scratch,*) this_field_p%i_value
3612 call strip_front_blanks(scratch)
3613 write (method(num_meth),'(a,a)') trim(method(num_meth)), &
3614 trim(this_field_p%name)
3615 write (control(num_meth),'(a)') &
3617 num_meth = num_meth + 1
3622 write (method(num_meth),'(a,a)') trim(method(num_meth)), &
3623 trim(this_field_p%name)
3624 write (control(num_meth),'(l1)') &
3625 this_field_p%l_value
3626 num_meth = num_meth + 1
3630 if(allocated(this_field_p%r_value)) write (scratch,*) this_field_p%r_value
3631 call strip_front_blanks(scratch)
3632 write (method(num_meth),'(a,a)') trim(method(num_meth)), &
3633 trim(this_field_p%name)
3634 write (control(num_meth),'(a)') &
3636 num_meth = num_meth + 1
3640 write (method(num_meth),'(a,a)') trim(method(num_meth)), &
3641 trim(this_field_p%name)
3642 write (control(num_meth),'(a)') &
3643 trim(this_field_p%s_value(1))
3644 do i = 2, this_field_p%max_index
3645 write (control(num_meth),'(a,a,$)') comma//trim(this_field_p%s_value(i))
3647 num_meth = num_meth + 1
3656 this_field_p => this_field_p%next
3660 end function find_method
3662 #include "field_manager_r4.fh"
3663 #include "field_manager_r8.fh"
3665 end module field_manager_mod
3667 ! close documentation grouping