FishMet: Fish feeding and appetite model, OPEN EDITION  0.1
FishMet: Fish feeding and appetite model, OPEN EDITION
m_common.f90 (16563)
Go to the documentation of this file.
1 !> @file m_common.f90
2 !! The FishMet Model: Definition of global parameters and
3 !! computational utilities.
4 
5 ! ------------------------------------------------------------------------------
6 ! Copyright notice:
7 ! Copyright (c) 2020-2024 Sergey Budaev & Ivar Rønnestad
8 ! FishMet model code is distributed according to GNU AGPL v3, see LICENSE.md
9 ! ------------------------------------------------------------------------------
10 
11 !> This module defines global parameters and general-level computational
12 !! utilities.
13 module commondata
14 
15  use realcalc
16  use, intrinsic :: iso_fortran_env, only : error_unit
17 
18  implicit none
19 
20  !! The include file is generated automatically from the global version
21  !! by the Make system. This functionality can be disabled by defining
22  !! an appropriate public string parameter `SVN_VERSION_GLOBAL`:
23  !! `character(len=*), parameter, public :: SVN_VERSION_GLOBAL = 'disabled'`
24  !!
25  !! Here is the Makefile macro that produces the version header
26  !! @verbatim
27  !! # Version string is obtained from SVN
28  !! VERSION_STRING := $(shell svn info --show-item last-changed-revision)
29  !! @endverbatim
30  !include 'version.h' ! autogenerated from svn macro in Makefile
31  character(len=*), parameter, public :: svn_version_global = '0.1_o'
32 
33  !> Program title that appears on the main interface/window.
34  character(len=*), parameter, public :: program_title = &
35  "FishMet: the fish appetite and feeding model"
36 
37  !> Logical flag that sets the DEBUG mode. When the program is running in
38  !! the **DEBUG mode**, additional checks are done and more output is
39  !! produced.
40  logical, parameter, public :: is_debug = .true.
41 
42  !> Compile-time logical flag that defines if stress effects described by
43  !! ::global_stress_factor_hour, ::global_stress_fact_suppress and
44  !! ::global_stress_cost_smr
45  logical, parameter, public :: is_stress_enable = .true.
46 
47  !> Warning level in the debugging mode.
48  real(srp), parameter, public :: debug_warn_level = 0.75_srp
49 
50  !> The default length of labels and similar text strings
51  integer, parameter, public :: label_len = 24
52 
53  !> The default length of the file name
54  integer, parameter, public :: file_name_len = 255
55 
56  !> The minimum length of a data file. Note that cancelled file name dialog
57  !! returns a six-char random sequence. Hence the minimum file name > 7
58  integer, parameter, public :: minimum_filename_chars = 7
59 
60  !> Global time constants, number of sec in hour and min
61  integer, parameter, public :: hour=3600, minute=60
62 
63  !> The maximum number of data points that are plotted for the model
64  !! output arrays (squeezed arrays).
65  integer, parameter, public :: max_points_plot = 1000
66 
67  !> Global parameter defining if progress bar is printed in the GUI mode.
68  logical, parameter, public :: is_progress_bar_gui = .true.
69 
70  !> Maximum value (steps) in the progress bar (minimum is 0)
71  integer, parameter, public :: gui_progress_bar_max = 100
72 
73  !> Number of steps in the progress bar at the GUI mode.
74  integer, parameter, public :: gui_progress_bar_steps = 100
75 
76  !> The symbol used for writing the progress bar in text mode.
77  character(len=1), parameter, public :: prog_bar_char = "#"
78 
79  !> The length of the progress bar in text mode.
80  integer, parameter, public :: prog_bar_len = 42
81 
82  ! New line character
83  character, parameter, public :: crlf = char(10)
84 
85  !> Global variable that defines the current time step.
86  integer(LONG), public :: global_time_step
87 
88  !> Default file base-name for saving model output arrays
89  character(len=*), parameter, public :: output_arrays_csv_file &
90  = "model_output_arrays"
91 
92  !> Default file base-name for saving model rate data
93  character(len=*), parameter, public :: output_rate_data_csv_file &
94  = "model_rate_data"
95 
96  !> Default output directory: all plot and data file will be saved into this
97  !! directory by default (local directory) unless reconfigured via
98  !! commondata::global_output_dest
99  character(len=*), parameter, public :: output_dest_dir = ""
100 
101  !> Output directory: all plot and output files will be saved into this
102  !! directory, default value is set by ::output_dest_dir.
103  !! Specifically, the following outputs are saved to the `output_dest`
104  !! directory:
105  !! - output plots produced by `save` command
106  !! (see if_cmd::cmd_iface_process_commands(), usually set by `setfil()`);
107  !! - output arrays by `save model_output`
108  !! (see simulate::output_arrays_save_csv());
109  !! - rate arrays by `save rate_output`
110  !! (see simulate::output_arrays_save_rate_data_csv());
111  !! - output statistics by `output stats` and `append stats`
112  !! (see if_cmd::cmd_iface_process_commands()).
113  !!
114  !! @note Note that fish object data `save fish_object_data` is saved into the
115  !! commondata::global_fish_object_binary_file exactly as specified in
116  !! the `fish_object_file` (commondata::global_fish_object_binary_file)
117  !! because this is non-human readable internal binary data not intended
118  !! as the model output.
119  !!
120  !! @note Note that this parameter affects the program behaviour only
121  !! in the command line user interface (if_cmd) and not in the graphical
122  !! user interface (if_gui).
123  !!
124  !! @note Note that the destination can refer to destination directory or
125  !! output file name prefix. In the former case it should terminate
126  !! with the operating system specific directory delimiter (e.g. `/` on
127  !! GNU/Linux or `/` on Microsoft Windows). If the value of
128  !! `output_dest` does not terminate this way, it will refer to the
129  !! file prefix.
130  !!
131  !! @verbatim
132  !! output_dest = "../../../SHARE/data/"
133  !! @endverbatim
134  character(len=FILE_NAME_LEN), public :: global_output_dest
135 
136  !> Defines the name of the file that keeps global parameters.
137  !! The parameter file contains the coupled pairs of the parameter
138  !! name and its value. Comments can be added, normally separated by
139  !! the `#` commet symbol
140  !! @verbatim
141  !! # Stomach mass of the fish, filling capacity
142  !! stomach_capacity = 123.45
143  !! # Water uptake, proportional to food item m
144  !! water_uptake = 0.1
145  !! @endverbatim
146  character(len=*), parameter, public :: param_file_name_def = "parameters.cfg"
147 
148  !> Defines the name of the file that keeps global parameters.
149  !! This name is normally defined by the environment variable
150  !! `FFA_MODEL_PARAMETER_FILE` that is processed in runtime::system_init().
151  !! If this environment variable is not set, the parameter file name is
152  !! defined by the default ::param_file_name_def.
153  character(len=FILE_NAME_LEN), public :: global_param_file_name
154 
155  ! Global parameters that are obtained from the parameter configuration file
156  ! Note that global parameters have GLOBAL_ prefix
157 
158  !> Maximum number of food items in the fish stomach, this is the maximum
159  !! index for the food items array in the_fish::stomach object.
160  !! @warning This index value must always be large enough to ensure maximum
161  !! filling capacity of the fish stomach ::global_stomach_mass in
162  !! units of food items ::Global_Food_Item_Mass given the mass
163  !! of many food items in the array is diminishing to zero.
164  !! MAX_FOOD_ITEMS_INDEX >> stomach_capacity / food_item_mass
165  !! @note Using a fixed size rather than allocatable array makes it easier
166  !! to program, faster runtime (less dynamic memory) and reduces
167  !! stack usage.
168  integer, parameter, public :: max_food_items_index = 4000
169 
170  !> Global logical flag to suppress extra text diagnostics and messages,
171  !! the "quiet mode" or the "verbose mode". If the parameter is TRUE, then
172  !! verbose mode is activated.
173  logical, public :: verbose
174 
175  !> Definition of time plot scale types for ::global_time_plots_unit_scale
176  enum, bind(C)
178  end enum
179 
180  !> @note Note that scale types are public.
182 
183  !> Default rate discretization interval in minutes, this interval is used
184  !! to plot the ingestion rate
185  !! @verbatim
186  !! rate_interval = 10
187  !! @endverbatim
188  integer, public :: global_rate_interval
189 
190  !> Default duration of the model run in **hours**. Note that parameter
191  !! file uses hours as unit. Conversion to time steps is done in
192  !! simulate::total_timesteps().
193  !! @verbatim
194  !! run_model_hours = 24
195  !! @endverbatim
196  integer, public :: global_run_model_hours
197 
198  !> Default hour at which the "daytime" is normally started. This means,
199  !! in particular, that at the start of the simulation day starts with the
200  !! offset equal to this parameter value.
201  !! @verbatim
202  !! day_starts_hour = 4
203  !! @endverbatim
204  integer, public :: global_day_starts_hour
205 
206  !> The offset (delay) to start feeding at the beginning of the simulation,
207  !! Note that the parameter file uses minutes as the unit.
208  !! @verbatim
209  !! feed_start_offset = 60
210  !! @endverbatim
211  integer, public :: global_run_model_feed_offset
212 
213  !> Default duration of the day time in hours, night duration is defined
214  !! as 24 - ::global_hours_day_feeding. Feeding activity occurs only
215  !! during the day time and not during the night.
216  !! @verbatim
217  !! daytime_hours = 12
218  !! @endverbatim
219  integer, public :: global_hours_daytime_feeding
220 
221  !> Fish body mass, g, at the start of the simulation
222  !! Configuration file example:
223  !! @verbatim
224  !! body_mass = 200.0
225  !! @endverbatim
226  real(srp), public :: global_body_mass
227 
228  !> Fish stomach mass, max. filling capacity.
229  !! Configuration file example:
230  !! @verbatim
231  !! stomach_capacity = 24.0
232  !! @endverbatim
233  real(srp), public :: global_stomach_mass
234 
235  !> Fish mid-gut mass capacity, max. filling capacity.
236  !! Configuration file example:
237  !! @verbatim
238  !! midgut_capacity = 24.0
239  !! @endverbatim
240  real(srp), public :: global_midgut_mass
241 
242  !> Logical flag that specifies that the stomach and midgut mass (capacity)
243  !! are automatically recalculated based on the body mass ::global_body_mass,
244  !! possibly overriding the parameter values set by ::global_stomach_mass and
245  !! ::global_midgut_mass
246  !! @verbatim
247  !! stomach_midgut_automatic = True
248  !! @endverbatim
250 
251  !> Maximum absorption ratio relative to the original dry food item
252  !! mass.
253  !! Configuration file example:
254  !! @verbatim
255  !! absorption_ratio = 0.5
256  !! @endverbatim
257  real(srp), public :: global_absorption_ratio
258 
259  !> Delay of ingestion, min
260  !! @note Note that the model is discrete, with a single time step
261  !! equivalent to one s
262  !! Configuration file example:
263  !! @verbatim
264  !! ingestion_delay = 30
265  !! @endverbatim
266  integer, public :: global_ingestion_delay_min
267 
268  !> Delay of digestion, min. It is the delay of absorption after a food item
269  !! was transmitted to the mid-gut. Any absorption can start only after
270  !! the digestion delay.
271  !! Configuration file example:
272  !! @verbatim
273  !! digestion_delay = 20
274  !! @endverbatim
275  integer, public :: global_digestion_delay_min
276 
277  !> The maximum duration a food item can be processed in the fish
278  !! mid-gut, min. If it stays in the mid-gut for any longer time, it is
279  !! excreted.
280  !! Configuration file example:
281  !! @verbatim
282  !! midgut_maxdur = 180
283  !! @endverbatim
285 
286  !> Proportion of water uptake, relative of dry mass of the food item.
287  !! Configuration file example:
288  !! @verbatim
289  !! water_uptake = 0.2
290  !! @endverbatim
291  real(srp), public :: global_water_uptake
292 
293  !> Parameters of the logistic water uptake function that defines the
294  !! temporary pattern of water uptake
295  !! `c /(1+a*%e**(-r2*t));` where c is the upper limit of the mass,
296  !! a and r are logistic parameters and %e is natural logarithm base.
297  !! Thus, water uptake dynamics is defined by two additional logistic
298  !! parameters: `Global_Water_Uptake_A` and `Global_Water_Uptake_R`
299  !! @note Maxima plot: `wxplot2d( 0.1/(1+200*%e**(-0.01*x)), [x,0,1900] );`
300  !! @verbatim
301  !! water_uptake_a = 200.0
302  !! @endverbatim
303  real(srp), public :: global_water_uptake_a
304 
305  !> Parameters of the logistic water uptake function that defines the
306  !! temporary pattern of water uptake
307  !! `c /(1+a*%e**(-r2*t));` where c is the upper limit of the mass,
308  !! a and r are logistic parameters and %e is natural logarithm base.
309  !! Thus, water uptake dynamics is defined by two additional logistic
310  !! parameters: `Global_Water_Uptake_A` and `Global_Water_Uptake_R`
311  !! @note Maxima plot: `wxplot2d( 0.1/(1+200*%e**(-0.01*x)), [x,0,1900] );`
312  !! @verbatim
313  !! water_uptake_r = 0.01
314  !! @endverbatim
315  real(srp), public :: global_water_uptake_r
316 
317  !> Food transport in stomach: Dimensionality of the transport pattern in the
318  !! stomach.
319  integer, public :: global_transport_dim
320  !> Food transport in stomach: Time grid array for interpolation.
321  !! Configuration file example:
322  !! @verbatim
323  !! transport_pattern_t = 0 2 8 12
324  !! @endverbatim
325  integer, allocatable, dimension(:), public :: global_transport_pattern_t
326  !> Food transport in stomach: Proportion of food mass left in stomach.
327  !! Configuration file example:
328  !! @verbatim
329  !! transport_pattern_r = 1.0 0.9 0.2 0.0
330  !! @endverbatim
331  real(srp), allocatable, dimension(:), public :: global_transport_pattern_r
332 
333  !> Temperature adjustment for absorption is controlled by the two
334  !! parameters below that define the temperature adfjustment grid factor
335  !! for the Michaelis-Menten equation: ::global_temp_factor_midgut_t and
336  !! ::global_temp_factor_midgut_m
337  !! This parameter provides the temperature grid for
338  !! Michaelis-Menten absorption process
339  !! @verbatim
340  !! midgut_temp_fact_t = c(5.0, 10.0, 15.0, 25.0)
341  !! @endverbatim
342  real(srp), allocatable, dimension(:), public :: global_temp_factor_midgut_t
343 
344  !> The pattern of stress effect on the appetite is described by two grid
345  !! arrays ::global_stress_factor_hour and
346  !! ::global_stress_fact_suppress
347  !! This ::global_stress_factor_hour defines the time grid for the
348  !! stress pattern, h
349  !! @verbatim
350  !! stress_grid_hour = c(0, 10, 20, 40, 50, 72)
351  !! @endverbatim
352  real(srp), allocatable, dimension(:), public :: global_stress_factor_hour
353 
354  !> The pattern of stress effect on the appetite is described by two grid
355  !! arrays ::global_stress_factor_hour and
356  !! ::global_stress_fact_suppress
357  !! This ::global_stress_factor_hour defines the time grid for the
358  !! stress pattern, h
359  !! @verbatim
360  !! stress_grid_fact = c(1.0, 0.9, 0.7, 0.2, 0.1, 0.0)
361  !! @endverbatim
362  real(srp), allocatable, dimension(:), public :: &
364 
365 
366  !> Maximum value of the metabolic cost of stress in units of resting metabolic
367  !! rate (SMR). The real metabilic cost scales depending on the time since
368  !! stress intervention following the stress pattern described by the
369  !! `stress_grid_hour` and `stress_grid_fact` grid.
370  !!
371  !! @verbatim
372  !! stress_metabolic_cost = 0.5
373  !! @endverbatim
374  real(srp), public :: global_stress_cost_smr
375 
376  !> Maximum value of the suppressive effect of stress on baseline activity.
377  !! The actual suppression of the baseline activity scales depending on the
378  !! time since stress intervention following the stress pattern described by
379  !! the `stress_grid_hour` and `stress_grid_fact` grid.
380  !! @verbatim
381  !! stress_inactivity = 0.5
382  !! @endverbatim
383  real(srp), public :: global_stress_activity_decr
384 
385  !> Global logical flag setting the time unit for the stress intervention
386  !! parameter array defined by ::global_stress_intervention_time. If set
387  !! to TRUE, the `stress` parameter is set in minutes, otherwise in raw
388  !! seconds.
389  !! @verbatim
390  !! stress_is_min = TRUE
391  !! @endverbatim
393 
394  !> The time of stress intervention is set by this array. It is the start
395  !! points of stressful effects (s).
396  !! @verbatim
397  !! # Two stress events, at 14h and 64h
398  !! # set in s
399  !! stress = c(50400, 240400)
400  !! # set in min
401  !! stress = c(840, 3840)
402  !! @endverbatim
403  integer(LONG), allocatable, dimension(:), public :: &
405 
406  !> Temperature adjustment for absorption is controlled by the two
407  !! parameters below that define the temperature adfjustment grid factor
408  !! for the Michaelis-Menten equation: ::global_temp_factor_midgut_t and
409  !! ::global_temp_factor_midgut_m
410  !! This parameter provides the temperature grid for
411  !! Michaelis-Menten absorption process
412  !! @verbatim
413  !! midgut_temp_fact_m = c(0.4, 0.6, 1.0, 8.0)
414  !! @endverbatim
415  real(srp), allocatable, dimension(:), public :: global_temp_factor_midgut_m
416 
417  !> The baseline *temperature* that applies to the stomach transport pattern
418  !! defined by ::global_transport_pattern_t and ::global_transport_pattern_r
419  !! @verbatim
420  !! transport_pattern_base_temp = 16
421  !! @endverbatim
423 
424  !> The baseline *fish mass* that applies to the stomach transport pattern
425  !! defined by ::global_transport_pattern_t and ::global_transport_pattern_r
426  !! @verbatim
427  !! transport_pattern_base_mass = 100
428  !! @endverbatim
430 
431  !> A data structure that defines the stomach emptying pattern as dependent
432  !! on the temperature and the fish mass
433  !! @verbatim
434  !! | 50 100 300 500 <--- grid_mass
435  !! -------+-------------------
436  !! 5 | 24 35 75 100
437  !! 10 | 15 20 50 75
438  !! 15 | 9 15 40 50
439  !! 20 | 5 10 30 35
440  !! 22 | 4 8 29 33
441  !! ^
442  !! +--- grid_temperature
443  !! @endverbatim
445  real(srp), dimension(:), allocatable :: grid_mass
446  real(srp), dimension(:), allocatable :: grid_temperature
447  integer, dimension(:,:), allocatable :: emptying_time
448  end type stomach_emptying_pattern
449 
450  !> Default stomach emptying pattern
452 
453  !> Fish appetite at night. This value is normally low because the fish do
454  !! not feed at night. If the value is absent or negative, fish appetite
455  !! does not reduce at night.
456  !! @verbatim
457  !! appetite_night=0.1
458  !! @endverbatim
459  real(srp), public :: global_appetite_fish_night
460 
461  !> Protective appetite threshold for stomach: this is the maximum value of
462  !! the the_fish::stomach::appetite_stomach() when the overall fish appetite
463  !! level the_fish::fish::appetite() depends only on stomach filling. Thereby
464  !! it can protect stomach from overfilling.
465  !! Configuration file example:
466  !! @verbatim
467  !! appetite_threshold_stomach = 0.2
468  !! @endverbatim
470 
471  !> Logistic function parameter **A** for the appetite factor.
472  !! See the_fish::appetite_func().
473  !! @verbatim
474  !! appetite_factor_a = 50000.0
475  !! @endverbatim
476  real(srp), public :: global_appetite_logist_a
477 
478  !> Logistic function parameter **R** for the appetite factor.
479  !! See the_fish::appetite_func().
480  !! @verbatim
481  !! appetite_factor_r = 20.0
482  !! @endverbatim
483  real(srp), public :: global_appetite_logist_r
484 
485  !> Michaelis-Meneten food absorption parameter in mid-gut, @f$ r_{max} @f$
486  !! (see the_fish::michaelis_menten()) relative to the maximum mass of the
487  !! food item, rate is per second (basic discrete step of the model).
488  !! @verbatim
489  !! midgut_michaelis_r_max = 0.8
490  !! @endverbatim
491  real(srp), public :: global_mid_gut_mm_r_max
492 
493  !> Michaelis-Meneten food absorption parameter in mid-gut, @f$ K_M @f$
494  !! (see the_fish::michaelis_menten()) relative to the total mid-gut
495  !! mass ::global_midgut_mass
496  !! @verbatim
497  !! midgut_michaelis_k = 0.25
498  !! @endverbatim
499  real(srp), public :: global_mid_gut_mm_k_m
500 
501  !> The steepness parameter of the Logistic energy component of appetite
502  !! @verbatim
503  !! appetite_energy_rate = 40.0
504  !! @endverbatim
505  real(srp), public :: global_energy_appetite_rate
506 
507  !> The shift parameter of the Logistic energy component of appetite
508  !! @verbatim
509  !! appetite_energy_shift = 0.2
510  !! @endverbatim
511  real(srp), public :: global_energy_appetite_shift
512 
513  !> Activity appetite factor determining how fish locomotor activity
514  !! increases with increasing appetite. See
515  !! the_fish::fish_activity_appetite_factor()
516  !! @verbatim
517  !! activity_appetite_factor = 0.5
518  !! @endverbatim
519  real(srp), public :: global_appetite_activity_factor
520 
521  !> A parameter defining branchial and urine (ZE+UE) energy consumption
522  !! as a factor to SMR, e.g. 0.3 means (UE+ZE) = 0.3*SMR
523  !! @note Note that this is the default method.
524  !! @verbatim
525  !! branchial_energy_factor = 0.3
526  !! @endverbatim
527  real(srp), public :: global_ue_ze_factor
528 
529  !> A parameter defining branchial and urine (ZE+UE) energy consumption
530  !! as a fixed ammonia excretion rate, micro (\mu) mol per g body mass
531  !! per hour
532  !! @verbatim
533  !! branchial_ammonia_rate = 0.5
534  !! @endverbatim
535  real(srp), public :: global_ue_ze_ammonia_excretion
536 
537  !> Logical flag to choose how branchial and urinal (ZE+UE) energy loss is
538  !! calculated.
539  !! - If this parameter is TRUE, fixed rate is assumed as defined by
540  !! ::global_ue_ze_ammonia_excretion
541  !! - if it is FALSE, a fraction of SMR defined by ::global_ue_ze_factor
542  !! is used for calculation. This method is the default.
543  !! .
544  logical, public :: global_is_ue_ze_fixed_rate
545 
546  !> Defines the specific dynamic action (SDA) that depends on the absorption
547  !! rate. There are two parameters that define the linear relationship between
548  !! absorption rate and SDA increase. This parameter defines the maximum
549  !! absorption rate when the maximum SDA (defined by `sda_energy_factor_max`)
550  !! is reached.
551  !! @verbatim
552  !! sda_absorption_rate_max = 0.0007
553  !! @endverbatim
554  real(srp), public :: global_sda_absorp_rate_max
555 
556  !> Defines the specific dynamic action (SDA) that depends on the absorption
557  !! rate. There are two parameters that define the linear relationship between
558  !! absorption rate and SDA increase. This parameter defines the maximum SDA
559  !! reached at the maximum absorption rate defined by
560  !! `sda_absorption_rate_max`.
561  !! For example, this parameter equal to 2.0 means that SDA = SMR * 2.0 at the
562  !! maximum absorption rate.
563  !! @verbatim
564  !! sda_energy_factor_max = 2.0
565  !! @endverbatim
566  real(srp), public :: global_sda_factor_max
567 
568  !> This parameter defines the X axis of the grid: temperature.
569  !! See commondata::oxygen_rate_std() for details. If the array is not
570  !! defined in the configuration file, use
571  !! commondata::trout_oxygen_grid_x_temp
572  !! Configuration file example:
573  !! @verbatim
574  !! smr_oxygen_temp = c(0.0, 0.5, 5.1, 15.2, 20.2, 25.2, 26.2)
575  !! @endverbatim
576  real(srp), allocatable, dimension(:), public :: global_oxygen_grid_x_temp
577 
578  !> This parameter defines the Y axis of the grid: oxygen consumption.
579  !! See commondata::oxygen_rate_std() for details. If the array is not
580  !! defined in the configuration file, use
581  !! commondata::trout_oxygen_grid_y_o2std
582  !! Configuration file example:
583  !! @verbatim
584  !! smr_oxygen_o2 = c(38.0, 38.0, 40.6, 67.8, 94.6, 136.7, 145.9)
585  !! @endverbatim
586  real(srp), allocatable, dimension(:), public :: global_oxygen_grid_y_o2std
587 
588  !> Standard mass of one food item.
589  !! Configuration file example:
590  !! @verbatim
591  !! food_item_mass = 2.0
592  !! @endverbatim
593  real(srp), public :: global_food_item_mass
594 
595  !> Gross energy content of the feed, MJ/kg (=kJ7g)
596  !! Configuration file example:
597  !! @verbatim
598  !! feed_gross_energy = 23.0
599  !! @endverbatim
600  real(srp), public :: global_food_gross_energy
601 
602  !> Minimum mass of food item that counts as non-zero,
603  !! a tolerance limit
604  real(srp), parameter, public :: min_mass_toler = 0.0001_srp
605 
606  !> Ambient temperature
607  !! Configuration file example:
608  !! @verbatim
609  !! temperature = 16.0
610  !! @endverbatim
611  real(srp), public :: global_temperature
612 
613  !> Standard rate of food item input, per minute.
614  !! Configuration file example:
615  !! @verbatim
616  !! food_input_rate = 6.0
617  !! @endverbatim
618  real(srp), public :: global_food_input_rate
619 
620  !> Parameters of the food provisioning pattern:
621  !! - interval food **is** provided, min
622  !! - interval food is **not** provided interval, min
623  !! .
624  !! These parameters determine the ::global_interval_food_pattern logical
625  !! array for each time step.
626  !! @verbatim
627  !! food_provision_pattern = 10 20
628  !! @endverbatim
629  integer, dimension(2), public :: global_interval_food_param
630 
631  !> Global food provisioning pattern, logical TRUE/FALSE for each time step.
632  logical, allocatable, dimension(:), public :: global_interval_food_pattern
633 
634  !> Global variable that keeps the file name for the food provisioning pattern
635  !! see environ::food_provisioning_get_file() for more details.
636  !! @note Note that the file has priority over the food provisioning pattern
637  !! defined in ::global_interval_food_param.
638  !! @verbatim
639  !! food_provision_file_name = "../zzz.csv"
640  !! @endverbatim
641  character(len=FILE_NAME_LEN), public :: global_food_pattern_file
642 
643  !> Global logical flag that defines if the feed scheduling pattern defined
644  !! by the ::global_interval_food_pattern file will propagate to all 24 h
645  !! periods or applied just once to the first such period.
646  !! @verbatim
647  !! food_provision_file_repeat = TRUE
648  !! @endverbatim
650 
651  !> Global logical flag that defines if the feed scheduling pattern defined
652  !! by the ::global_food_pattern_file file represents the data for each time
653  !! step (s) (if TRUE) or by minute (FALSE). The default value is FALSE i.e.
654  !! data are given for each minute.
655  !! @verbatim
656  !! food_provisioning_file_by_s = FALSE
657  !! @endverbatim
659 
660  !> Global baseline locomotor activity pattern (swimming speed) for each
661  !! time step consists of two components:
662  !! - baseline that differs at day ::global_baseline_activity_day and
663  !! night ::global_baseline_activity_night
664  !! - feeding activity increment defined by feeding activity factor
665  !! ::global_feeding_activity_factor
666  !! .
667  !!
668  !! This is illustrated by the following diagram:
669  !!
670  !! @verbatim
671  !! +----+ feeding activity multiplier
672  !! | |
673  !! +--+ - -+----------------+
674  !! | | baseline
675  !! -------------------+ +-------------------> activity
676  !! -night-------------+-day--------------------+-night------------->
677  !! @endverbatim
678  !!
679  !! see the_fish::fish_activity_baseline() for more details.
680  !!
681  !! @verbatim
682  !! baseline_activity_day = 4.0
683  !! baseline_activity_night = 0.5
684  !! @endverbatim
685  real(srp), public :: global_baseline_activity_day, &
687 
688  !> Locomotor activity multiplier factor for the activity during the feeding
689  !! time, i.e. when ::global_interval_food_pattern is TRUE.
690  real(srp), public :: global_feeding_activity_factor
691 
692  !> Global variable keeping the file name for saving general simulation
693  !! statistics. Data (rows) for each simulation run are normally appended
694  !! to this file. The data format is CSV.
695  !! See simulate::model_output_stats_row_csv() for some details.
696  !! @verbatim
697  !! stats_output_file = fishmet_stats.csv
698  !! @endverbatim
699  character(len=FILE_NAME_LEN), public :: global_stats_output_file
700 
701 
702  !> Global variable keeping the file name for the baseline stomach emptying
703  !! matrix file that keeps the stomach emptying times for fish of
704  !! different mass and at different temperatures. The file should be
705  !! prepared in the CSV format and have the following structure:
706  !! TODO
707  !! @verbatim
708  !! stomach_emptying_matrix = fishmet_stom_data.csv
709  !! @endverbatim
710  character(len=FILE_NAME_LEN), public :: global_stomach_emptying_matrix_file
711 
712  !> Define if the output stats are saved using the long or short format. In
713  !! the first case, both the input parameters and the output stats are
714  !! saved/appended into the CSV file. Note that the default value is TRUE.
715  !! @verbatim
716  !! stats_output_long = TRUE
717  !! @endverbatim
718  logical, public :: global_output_stats_is_long
719 
720  !> Default value of the general simulation statistics file defined in
721  !! commondata::global_stats_output_file
722  character(len=*), parameter, public :: def_stats_output_file = &
723  "fishmet_stats.csv"
724  !> Natural logarithm base, `e`
725  real(srp), parameter, public :: base_e = 2.7182818284590451_srp
726 
727  ! *Rainbow trout* default parameters that are normally unchanged in the model
728 
729  !> Interpolation grid defining the function calculating standard oxygen
730  !! consumption in rainbow trout based on the curve from Fig. 9 in
731  !! Evans, D.O. (1990) Metabolic Thermal compensation by rainbow trout:
732  !! effects on standard metabolic rate and potential usable power.
733  !! Trans. Am. Fish. Soc. 119, 585–600-
734  !!
735  !! Interpolation data from Evan 1990:
736  !!@verbatim
737  !! [0.0 0.5 5.1 15.2 20.2 25.2 26.2] \
738  !! [38.0 38.0 40.6 67.8 94.6 136.7 145.9]
739  !!@endverbatim
740  !!
741  !! This parameter defines the X axis of the grid: temperature.
742  !! See commondata::oxygen_rate_std() for details.
743  real(srp), dimension(*), parameter :: trout_oxygen_grid_x_temp = &
744  [ 0.0_srp, 0.5_srp, 5.0_srp, 15.2_srp, 20.1_srp, 25.2_srp, 26.2_srp]
745 
746  !> Interpolation grid defining the function calculating standard oxygen
747  !! consumption in rainbow trout based on the curve from Fig. 9 in
748  !! Evans, D.O. (1990) Metabolic Thermal compensation by rainbow trout:
749  !! effects on standard metabolic rate and potential usable power.
750  !! Trans. Am. Fish. Soc. 119, 585–600-
751  !!
752  !! Interpolation data from Evan 1990:
753  !!@verbatim
754  !! [0.0 0.5 5.1 15.2 20.2 25.2 26.2] \
755  !! [38.0 38.0 40.6 67.8 94.6 136.7 145.9]
756  !!@endverbatim
757  !!
758  !! This parameter defines the Y axis of the grid: oxygen consumption.
759  !! See commondata::oxygen_rate_std() for details.
760  real(srp), dimension(*), parameter :: trout_oxygen_grid_y_o2std = &
761  [38.0_srp, 38.0_srp, 40.6_srp, 67.8_srp, 94.6_srp, 136.7_srp, 145.9_srp]
762 
763  !> Time scale abbreviation strings for Y axes of output plots
764  character(len=*), parameter, public :: txt_sec="s",txt_min="min",txt_hour="h"
765 
766  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
767 
768  !> Public interface for the cum2rate function that calculates rate from
769  !! raw cumulative array. Note that rate unit time is defined by the
770  !! `integer` rate_unit argument (as default, per commondata::minute).
771  !! For details see ::cum2rate_int() and ::cum2rate_r().
772  interface cum2rate
773  module procedure cum2rate_int
774  module procedure cum2rate_r
775  end interface cum2rate
776 
777  !-----------------------------------------------------------------------------
778  !> Simple history stack function, add to the end of the stack. We need
779  !! only to add components on top (end) of the stack and retain
780  !! `HISTORY_SIZE_SPATIAL` elements of the prior history (for a spatial moving
781  !! object). The stack works as follows, assuming 100 and 200 are added:
782  !! @n [1 2 3 4 5 6 7 8 9 10]
783  !! @n [2 3 4 5 6 7 8 9 10 **100**]
784  !! @n [3 4 5 6 7 8 9 10 100 **200**]
785  interface add_to_history
786  module procedure add_to_history_i4
787  module procedure add_to_history_r
788  module procedure add_to_history_char
789  end interface add_to_history
790 
791  !> Return the last element of a history array with optional offset.
792  interface last
793  module procedure last_in_history_i4
794  module procedure last_in_history_r
795  module procedure last_in_history_char
796  end interface last
797 
798  !> Function calculating standard oxygen consumption in rainbow trout
799  !! based on the curve from Fig. 9 in
800  !! Evans, D.O. (1990) Metabolic Thermal compensation by rainbow trout:
801  !! effects on standard metabolic rate and potential usable power.
802  !! Trans. Am. Fish. Soc. 119, 585–600
803  !! @note Note that `oxygen_rate_std` is implemented in both array
804  !! and scalar versions.
805  interface oxygen_rate_std
806  module procedure oxygen_rate_std_spline
807  module procedure oxygen_rate_std_ddpi
808  end interface oxygen_rate_std
809 
810  !> Function for calculating specific growth rate
811  interface sgr
812  module procedure sgr_r
813  module procedure sgr_i
814  module procedure sgr_blockrate
815  end interface
816 
817  !> Force a value within the range set by the vmin and vmax dummy parameter
818  !! values.
819  interface within
820  module procedure within_i
821  module procedure within_r
822  end interface within
823  ! @note Note that all component procedures are private.
824 
825  ! @note Note that all component procedures are private.
830 
831 contains
832 
833  !> Produce a text string containing all global parameters of the model
834  impure function model_parameters_txt () result ( message_show )
835  use base_utils, only: tostr
836 
837  character(len=:), allocatable :: message_show
838 
839  character(len=500) :: msg_line
840  character(len=LABEL_LEN) :: yes_no_descr, f_min_s_descr, f_exist_descr
841  logical :: is_okay
842 
843  character(len=*), parameter :: fmt_real="f9.3",fmt_rsrt="f8.2",fmt_int="I9"
844  character(len=*), parameter :: fmt_real_lon="f14.3"
845  character(len=:), allocatable :: fmtstr, label_stats
846 
847  message_show = crlf // " GLOBAL PARAMETERS" // crlf
848 
849  fmtstr = "(a,t56," // fmt_int // ")"
850  write(msg_line,fmtstr) " Run duration, hours ('run_model_hours'): ", &
852  message_show = message_show // crlf // trim(msg_line)
853 
854  fmtstr = "(a,t56," // fmt_int // ")"
855  write(msg_line,fmtstr) " Day start offset, hours ('day_starts_hour'): ", &
857  message_show = message_show // crlf // trim(msg_line)
858 
859  fmtstr = "(a,t56," // fmt_int // ")"
860  write(msg_line,fmtstr) " Daytime duration, hours ('daytime_hours'): ", &
862  message_show = message_show // crlf // trim(msg_line)
863 
864  fmtstr = "(a,t56," // fmt_real // ")"
865  write(msg_line,fmtstr) " Temperature ('temperature'): ", &
867  message_show = message_show // crlf // trim(msg_line)
868 
869  fmtstr = "(a,t56," // fmt_real // ")"
870  write(msg_line,fmtstr) " Body mass ('body_mass'): ", global_body_mass
871  message_show = message_show // crlf // trim(msg_line)
872 
873  fmtstr = "(a,t56," // fmt_real // ")"
874  write(msg_line,fmtstr) " Stomach mass ('stomach_capacity'): ", &
876  message_show = message_show // crlf // trim(msg_line)
877 
878  fmtstr = "(a,t56," // fmt_real // ")"
879  write(msg_line,fmtstr) " Mid-gut mass ('midgut_capacity'): ", &
881  message_show = message_show // crlf // trim(msg_line)
882 
884  label_stats = " YES"
885  else
886  label_stats = " NO"
887  end if
888  message_show = message_show // crlf // &
889  " Stomach/midgut capacity auto (`stomach_midgut_automatic`): " // &
890  label_stats
891 
892  fmtstr = "(a,t56," // fmt_real // ")"
893  write(msg_line,fmtstr) " Absorption ratio ('absorption_ratio'): ", &
895  message_show = message_show // crlf // trim(msg_line)
896 
897  fmtstr = "(a,t56," // fmt_int // ")"
898  write(msg_line,fmtstr) " Delay after ingestion, min (ingestion_delay): ", &
900  message_show = message_show // crlf // trim(msg_line)
901 
902  fmtstr = "(a,t56," // fmt_real // ")"
903  write(msg_line,fmtstr) " Water Uptake ('water_uptake'): ", &
905  message_show = message_show // crlf // trim(msg_line)
906 
907  fmtstr = "(a,t56," // fmt_real // ")"
908  write(msg_line,fmtstr) " Water Uptake logistic a ('water_uptake_a'): ", &
910  message_show = message_show // crlf // trim(msg_line)
911 
912  fmtstr = "(a,t56," // fmt_real // ")"
913  write(msg_line,fmtstr) " Water Uptake logistic r ('water_uptake_r'): ", &
915  message_show = message_show // crlf // trim(msg_line)
916 
917  write(msg_line,"(a)") " Transport pattern in stomach " // &
918  "('transport_pattern_t','transport_pattern_r'):"
919  message_show = message_show // crlf // trim(msg_line)
920 
921  fmtstr="(a,t14,a,"//tostr(size(global_transport_pattern_t))//"I8"//",a)"
922  write(msg_line,fmtstr) " Grid X: ", "[", global_transport_pattern_t, " ]"
923  message_show = message_show // crlf // trim(msg_line)
924 
925  fmtstr="(a,t14,a,"//tostr(size(global_transport_pattern_r))//fmt_rsrt//",a)"
926  write(msg_line,fmtstr) " Grid Y: ", "[", global_transport_pattern_r, " ]"
927  message_show = message_show // crlf // trim(msg_line)
928 
929  fmtstr = "(a,t51," // fmt_real_lon // ")"
930  write(msg_line,fmtstr) " Appetite logistic a ('appetite_factor_a'): ", &
932  message_show = message_show // crlf // trim(msg_line)
933 
934  fmtstr = "(a,t56," // fmt_real // ")"
935  write(msg_line,fmtstr) " Appetite logistic r ('appetite_factor_r'): ", &
937  message_show = message_show // crlf // trim(msg_line)
938 
939  fmtstr = "(a,t56," // fmt_real // ")"
940  write(msg_line,fmtstr) &
941  " Appetite at night ('appetite_night'): ", &
943  message_show = message_show // crlf // trim(msg_line)
944 
945  fmtstr = "(a,t56," // fmt_real // ")"
946  write(msg_line,fmtstr) &
947  " Appetite threshold ('appetite_threshold_stomach'): ", &
949  message_show = message_show // crlf // trim(msg_line)
950 
951  fmtstr = "(a,t56," // fmt_real // ")"
952  write(msg_line,fmtstr) &
953  " Activity appetite factor ('activity_appetite_factor'): ",&
955  message_show = message_show // crlf // trim(msg_line)
956 
957  fmtstr = "(a,t56," // fmt_int // ")"
958  write(msg_line,fmtstr) " Digestion Delay, min ('digestion_delay'): ", &
960  message_show = message_show // crlf // trim(msg_line)
961 
962  fmtstr = "(a,t56," // fmt_int // ")"
963  write(msg_line,fmtstr) " Max. duration in mid-gut, h ('midgut_maxdur'): ",&
965  message_show = message_show // crlf // trim(msg_line)
966 
967  fmtstr = "(a,t56," // fmt_real // ")"
968  write(msg_line,fmtstr) &
969  " Energy appetite steepness ('appetite_energy_rate'): ", &
971  message_show = message_show // crlf // trim(msg_line)
972 
973  fmtstr = "(a,t56," // fmt_real // ")"
974  write(msg_line,fmtstr) &
975  " Energy appetite steepness ('appetite_energy_shift'): ", &
977  message_show = message_show // crlf // trim(msg_line)
978 
979  write(msg_line,"(a)") " Mid-gut Michaelis-Meneten absorption constants:"
980  message_show = message_show // crlf // trim(msg_line)
981 
982  fmtstr = "(a,t56," // "f9.7" // ")"
983  write(msg_line,fmtstr) " r_max ('midgut_michaelis_r_max'): ", &
985  message_show = message_show // crlf // trim(msg_line)
986 
987  fmtstr = "(a,t56," // fmt_real // ")"
988  write(msg_line,fmtstr) " K_M ('midgut_michaelis_k'): ", &
990  message_show = message_show // crlf // trim(msg_line)
991 
992  write(msg_line,"(a)") " Michaelis-Meneten temperature adjustment grid:"
993  message_show = message_show // crlf // trim(msg_line)
994 
995  fmtstr = "(a,t26,a," // tostr(size(global_temp_factor_midgut_t)) // &
996  fmt_rsrt // ",a)"
997  write(msg_line,fmtstr) " 'midgut_temp_fact_t': ", "[", &
999  message_show = message_show // crlf // trim(msg_line)
1000  !
1001  fmtstr = "(a,t26,a," // tostr(size(global_temp_factor_midgut_m)) // &
1002  fmt_rsrt // ",a)"
1003  write(msg_line,fmtstr) " 'midgut_temp_fact_m': ", "[", &
1005  message_show = message_show // crlf // trim(msg_line)
1006 
1007  write(msg_line,"(a)") " SMR function ('smr_oxygen_temp', 'smr_oxygen_o2'):"
1008  message_show = message_show // crlf // trim(msg_line)
1009 
1010  fmtstr = "(a,t14,a," // tostr(size(global_oxygen_grid_x_temp)) // &
1011  fmt_rsrt // ",a)"
1012  write(msg_line,fmtstr) " Grid X: ", "[", global_oxygen_grid_x_temp, " ]"
1013  message_show = message_show // crlf // trim(msg_line)
1014  !
1015  fmtstr = "(a,t14,a," // tostr(size(global_oxygen_grid_y_o2std)) // &
1016  fmt_rsrt // ",a)"
1017  write(msg_line,fmtstr) " Grid Y: ", "[", global_oxygen_grid_y_o2std, " ]"
1018  message_show = message_show // crlf // trim(msg_line)
1019 
1020  if (global_is_ue_ze_fixed_rate) then
1021  fmtstr = "(a,t56," // fmt_real // ")"
1022  write(msg_line,fmtstr) &
1023  " Branchial/urine ('branchial_ammonia_rate'): ", &
1025  message_show = message_show // crlf // trim(msg_line)
1026  else
1027  fmtstr = "(a,t56," // fmt_real // ")"
1028  write(msg_line,fmtstr) &
1029  " Branchial/urine ('branchial_energy_factor'): ", &
1031  message_show = message_show // crlf // trim(msg_line)
1032  end if
1033 
1034  fmtstr = "(a,t57," // "f8.6" // ")"
1035  write(msg_line,fmtstr) &
1036  " SDA absorp. ('sda_absorption_rate_max'): ", &
1038  message_show = message_show // crlf // trim(msg_line)
1039 
1040  fmtstr = "(a,t56," // fmt_real // ")"
1041  write(msg_line,fmtstr) &
1042  " SDA max level ('sda_energy_factor_max'): ", &
1044  message_show = message_show // crlf // trim(msg_line)
1045 
1046  fmtstr = "(a,t56," // fmt_real // ")"
1047  write(msg_line,fmtstr) &
1048  " Baseline day activity ('baseline_activity_day'): ", &
1050  message_show = message_show // crlf // trim(msg_line)
1051 
1052  fmtstr = "(a,t56," // fmt_real // ")"
1053  write(msg_line,fmtstr) &
1054  " Baseline night activity ('baseline_activity_night'): ",&
1056  message_show = message_show // crlf // trim(msg_line)
1057 
1058  fmtstr = "(a,t56," // "f9.4" // ")"
1059  write(msg_line,fmtstr) " Food item mass ('food_item_mass'): ", &
1061  message_show = message_show // crlf // trim(msg_line)
1062 
1063  fmtstr = "(a,t40," // fmt_real // ",a," // "I3" // ",a)"
1064  write(msg_line,fmtstr) " Food input rate ('food_input_rate'): ", &
1066  "/min = every ", rate2int(global_food_input_rate), " s"
1067  message_show = message_show // crlf // trim(msg_line)
1068 
1069  ! Print either food provision pattern or provision data file
1070  if ( global_food_pattern_file == "" ) then
1071 
1072  fmtstr = "(a,t56,a," // tostr(size(global_interval_food_param)) // "i8" &
1073  // ",a)"
1074  write(msg_line,fmtstr) " Food scheduling ('food_provision_pattern'): ", &
1075  "[", global_interval_food_param, " ]"
1076  message_show = message_show // crlf // trim(msg_line)
1077 
1078  fmtstr = "(a,t56," // fmt_int // ")"
1079  write(msg_line,fmtstr) " Offset to start feeding, min " // &
1080  "('feed_start_offset'):", global_run_model_feed_offset
1081  message_show = message_show // crlf // trim(msg_line)
1082 
1083  else
1084 
1086  yes_no_descr = " all 24h periods)"
1087  else
1088  yes_no_descr = " apply once)"
1089  end if
1090 
1092  f_min_s_descr = " (by s"
1093  else
1094  f_min_s_descr = " (by min"
1095  end if
1096 
1097  inquire( file=global_food_pattern_file, exist=is_okay )
1098 
1099  if ( is_okay ) then
1100  f_exist_descr = ""
1101  else
1102  f_exist_descr = " NO FILE"
1103  end if
1104 
1105  write(msg_line,"(5a)") " Food provision scheduling from file: ", &
1106  trim(global_food_pattern_file), trim(f_min_s_descr), &
1107  trim(yes_no_descr), trim(f_exist_descr)
1108  message_show = message_show // crlf // trim(msg_line)
1109 
1110  if ( .not. is_okay ) then
1111  fmtstr = "(a,t30,a," // tostr(size(global_interval_food_param)) // &
1112  "i8" // ",a)"
1113  write(msg_line,fmtstr) " Food provisioning schedule: ", "[", &
1115  message_show = message_show // crlf // trim(msg_line)
1116 
1117  fmtstr = "(a,t30," // fmt_int // ")"
1118  write(msg_line,fmtstr) " Offset to start feeding, min:", &
1120  message_show = message_show // crlf // trim(msg_line)
1121  end if
1122 
1123  end if
1124 
1125  if ( is_stress() ) then
1126  write(msg_line,"(a)") " Stress parameters:"
1127  message_show = message_show // crlf // trim(msg_line)
1128 
1129  fmtstr = "(a,t56," // fmt_real // ")"
1130  write(msg_line,fmtstr) " 'stress_metabolic_cost': ", &
1132  message_show = message_show // crlf // trim(msg_line)
1133 
1134  fmtstr = "(a,t56," // fmt_real // ")"
1135  write(msg_line,fmtstr) " 'stress_inactivity': ", &
1137  message_show = message_show // crlf // trim(msg_line)
1138 
1139  fmtstr = "(a,t26,a," // tostr(size(global_stress_intervention_time)) // &
1140  fmt_int // ",a)"
1142  write(msg_line,fmtstr) " 'stress' (min): ", "[", &
1144  else
1145  write(msg_line,fmtstr) " 'stress' (s): ", "[", &
1147  end if
1148  message_show = message_show // crlf // trim(msg_line)
1149 
1150  fmtstr = "(a,t26,a," // tostr(size(global_stress_factor_hour)) // &
1151  fmt_real // ",a)"
1152  write(msg_line,fmtstr) " 'stress_grid_hour: ", "[", &
1154  message_show = message_show // crlf // trim(msg_line)
1155 
1156  fmtstr = "(a,t26,a," // tostr(size(global_stress_fact_suppress)) // &
1157  fmt_real // ",a)"
1158  write(msg_line,fmtstr) " 'stress_grid_fact: ", "[", &
1160  message_show = message_show // crlf // trim(msg_line)
1161  end if
1162 
1163  message_show = message_show // crlf
1164 
1165  message_show = message_show // crlf // " Config file: '" // &
1166  trim(global_param_file_name) // "'"
1167 
1168  message_show = message_show // crlf // &
1169  " Output destination ('output_dest'): '" // &
1170  trim(global_output_dest) // "'"
1171 
1172  if (global_output_stats_is_long) then
1173  label_stats = " [long]"
1174  else
1175  label_stats = " [short]"
1176  end if
1177 
1178  message_show = message_show // crlf // " Global output stats file: '" // &
1179  trim(global_stats_output_file) // "'" // label_stats
1180 
1181  if (is_debug) write(msg_line,"(a)") " Running in the DEBUG mode;"
1182  message_show = message_show // crlf // crlf // trim(msg_line)
1183 
1184  fmtstr = "(2a)"
1185  write(msg_line,fmtstr) " Version information: ", svn_version_global
1186  message_show = message_show // trim(msg_line)
1187 
1188  message_show = message_show // crlf
1189 
1190  end function model_parameters_txt
1191 
1192  !> Calculate if stress factor grid is enabled. Stress is disabled if
1193  !! *any* of the stress parameters are *not defined*; this includes the
1194  !! stress suppression grid array or the stress intervention array
1195  pure function is_stress() result (is_enabled)
1196  logical :: is_enabled
1197 
1198  if ( (.not. is_stress_enable) .or. &
1199  any(global_stress_factor_hour == missing) .or. &
1200  any(global_stress_fact_suppress == missing) .or. &
1201  any(global_stress_intervention_time == unknown) .or. &
1202  global_stress_cost_smr == missing .or. &
1203  global_stress_activity_decr == missing ) then
1204  is_enabled = .false.
1205  else
1206  is_enabled = .true.
1207  end if
1208 
1209  end function is_stress
1210 
1211  !> Logistic function is defined by the equation
1212  !! @f[ \frac{c}{1+a \cdot e^{-r \cdot x}} @f]
1213  elemental function logistic( X, A, C, R ) result (Y)
1214  real(srp), intent(in) :: x
1215  real(srp), intent(in) :: a
1216  real(srp), intent(in) :: c
1217  real(srp), intent(in) :: r
1218  real(srp) :: y
1219 
1220  y = c / ( 1.0_srp + a * base_e **(-r*x) )
1221 
1222  end function logistic
1223 
1224  !> Michaelis-Meneten function:
1225  !! @f[ y = \frac{r_{max} x}{K_M + x} @f]
1226  !! Here @f$ r_{max} @f$ is the maximum rate, @f$ K_M @f$ is the Michaelis
1227  !! constant equal to the x at which 1/2 of the maximum rate @f$ r_{max} @f$
1228  !! is reached.
1229  !!
1230  !! *gnupolot code:*
1231  !! @verbatim
1232  !! set xrange [0:1]
1233  !! plot 0.8*x/(0.25+x)
1234  !! @endverbatim
1235  elemental function michaelis_menten( x, r_max, x50 ) result ( r )
1236  !> Independent variable.
1237  real(srp), intent(in) :: x
1238  !> @f$ r_{max} @f$ constant.
1239  real(srp), intent(in) :: r_max
1240  !> The Michaelis constant @f$ K_M @f$.
1241  real(srp), intent(in) :: x50
1242  !> Output response rate.
1243  real(srp) :: r
1244 
1245  r = r_max * x / (x50 + x)
1246 
1247  end function michaelis_menten
1248 
1249  !> Calculate temperature adjustment factor for the Michaelis-Menten
1250  !! absorption process in the midgut.
1251  !! The function is based on the cubic spline approcimation using the
1252  !! grid vectors defined by commondata::global_temp_factor_midgut_t and
1253  !! commondata::global_temp_factor_midgut_t.
1254  elemental function midgut_temp_factor(temperature) result (afact)
1255  use base_utils, only : cspline_scalar
1256  real(srp), optional, intent(in) :: temperature
1257  real(srp) :: afact
1258 
1259  real(srp) :: temp_loc
1260 
1261  if (present(temperature)) then
1262  temp_loc = temperature
1263  else
1264  temp_loc = global_temperature
1265  end if
1266 
1267  afact = cspline_scalar( global_temp_factor_midgut_t, &
1268  global_temp_factor_midgut_m, temp_loc )
1269 
1270  end function midgut_temp_factor
1271 
1272  !> Calculate the energy content of a feed mass.
1273  !! @note Units correspondence: mass in kg corresponds to energy in MJ,
1274  !! mass in g corresponds to kJ.
1275  elemental function feed_energy (feed_mass, gross_energy, &
1276  is_adjust_water_uptake, water_uptake) result (energy)
1277  !> The mass of the feed, kg
1278  real(srp), intent(in) :: feed_mass
1279  !> Gross energy content of the feed, MJ/kg (=kJ/g)
1280  real(srp), intent(in), optional :: gross_energy
1281  !> Optional logical parameter indicating if the energy content of the
1282  !! feed is adjusted for water uptake. If TRUE, the energy content of feed
1283  !! @f$ F( \Delta a ) @f$ is adjusted by subtracking the water added:
1284  !! @f$ F( \Delta a ) = E_G \Delta a - E_G \Delta a \times u @f$.
1285  !! The default value is TRUE.
1286  logical, intent(in), optional :: is_adjust_water_uptake
1287  !> Optional water uptake coefficient, if absent, the degault value is
1288  !! taken from ::global_water_uptake.
1289  real(srp), intent(in), optional :: water_uptake
1290  !> Energy content of the input feed mass, MJ
1291  real(srp) :: energy
1292 
1293  logical:: is_adjust_loc
1294  real(srp) :: water_uptake_loc
1295 
1296  if (present(is_adjust_water_uptake)) then
1297  is_adjust_loc = is_adjust_water_uptake
1298  else
1299  is_adjust_loc = .true.
1300  end if
1301 
1302  if (present(water_uptake)) then
1303  water_uptake_loc = water_uptake
1304  else
1305  water_uptake_loc = global_water_uptake
1306  end if
1307 
1308  ! @f$ \Delta a ) = E_G \Delta a @f$
1309  if (present(gross_energy)) then
1310  energy = feed_mass * gross_energy
1311  else
1312  energy = feed_mass * global_food_gross_energy
1313  end if
1314 
1315  ! If adjustment for uptake is used,
1316  ! @f$ \Delta a ) = E_G \Delta a - E_G \Delta a \times u @f$
1317  if (is_adjust_loc) then
1318  energy = energy - energy * water_uptake_loc
1319  end if
1320 
1321  end function feed_energy
1322 
1323  !> Calculate the enegretic equivalent (kJ) of Oxygen uptake (l O2)
1324  !! Source: Rønnestad B291.
1325  elemental function energy_o2(l_o2) result (energy_kj)
1326  !> Mass of oxygen, l
1327  real(srp), intent(in) :: l_o2
1328  !> Energy equivalent, kJ
1329  real(srp) :: energy_kj
1330 
1331  energy_kj = 19.4_srp * l_o2
1332 
1333  end function energy_o2
1334 
1335  !> Calculate the mass of 1 dm3 (l) of Oxygen.
1336  !!
1337  !! Based on this: 1 mol (22.4l) O2 is 32g, therefore, 1l = 1.43g
1338  elemental function o2mass(vol_l) result (mass_g)
1339  !> Volume of oxygen, l
1340  real(srp), intent(in) :: vol_l
1341  !> Mass of oxygen, g
1342  real(srp) :: mass_g
1343  real(srp), parameter :: ratio = 32.0_srp / 22.4_srp ! 1.4285714
1344 
1345  mass_g = ratio * vol_l
1346 
1347  end function o2mass
1348 
1349  !> Calculate the volume of 1 dm3 (l) of oxygen.
1350  elemental function o2vol(mass_g) result (vol_l)
1351  !> Mass of oxygen, g
1352  real(srp), intent(in) :: mass_g
1353  !> Volume of oxygen, l
1354  real(srp) :: vol_l
1355  real(srp), parameter :: ratio = 22.4_srp / 32.0_srp ! 0.7
1356 
1357  vol_l = mass_g * ratio
1358 
1359  end function o2vol
1360 
1361  !> The function converts rate per min to interval in seconds.
1362  elemental function rate2int(rate_min) result (n_sec)
1363  !> Optional rate per minute. Default value is set by
1364  !! commondata::global_food_input_rate.
1365  real(srp), optional, intent(in) :: rate_min
1366  !> returns interval in seconds between successive events (i.e. every N sec).
1367  integer(LONG) :: n_sec
1368 
1369  real(srp) :: rate_min_loc
1370 
1371  if (present(rate_min)) then
1372  rate_min_loc = rate_min
1373  else
1374  rate_min_loc = global_food_input_rate
1375  end if
1376 
1377  !> The main equation is @f$ n = \frac{60}{r_{min}} @f$, where
1378  !! @f$ n @f$ is the interval (s) and @f$ r_{min} @f$ is the rate
1379  !! per minute.
1380  !!
1381  !! Then, to do something every n seconds is implemented as
1382  !! @code
1383  !! if ( mod( time_step, rate2int( rate ) )==0 ) then
1384  !! @endcode
1385  n_sec = nint( minute/rate_min_loc, long )
1386 
1387  end function rate2int
1388 
1389 
1390  !> Given an in put array of cumulative sum, this function calculates the
1391  !! rate per unit time defined by the rate_unit argument (as default,
1392  !! per commondata::minute). This requires working on the appropriate
1393  !! difference values for each time interval.
1394  !! Note: This function version accepts *integer* array.
1395  pure function cum2rate_int( inarray, rate_unit ) result(rate_out)
1396  !> Input Output_Arrays array (integer type)
1397  integer, dimension(:), intent(in) :: inarray
1398  !> Optional rate parameter, the rate is calculated per this number of
1399  !! seconds, the default is per minute
1400  integer, intent(in), optional :: rate_unit
1401 
1402  real(srp), dimension(:), allocatable :: rate_out
1403 
1404  integer :: rate_unit_loc, oarray_max, sum_rate, i, j
1405 
1406  if (present(rate_unit)) then
1407  rate_unit_loc = rate_unit
1408  else
1409  rate_unit_loc = minute
1410  end if
1411 
1412  oarray_max = size(inarray) ! Output_Arrays
1413  allocate( rate_out(oarray_max/rate_unit_loc) )
1414 
1415  j = 0
1416  sum_rate = inarray(1)
1417 
1418  do i=2, oarray_max
1419  if( inarray(i) > inarray(i-1) ) then
1420  sum_rate = sum_rate + inarray(i) - inarray(i-1)
1421  end if
1422  if (mod(i,rate_unit_loc)==0) then
1423  j = j + 1
1424  rate_out(j) = real(sum_rate,srp) / real(rate_unit_loc,srp)
1425  sum_rate = 0
1426  end if
1427  end do
1428 
1429  end function cum2rate_int
1430 
1431  !> Given an in put array of cumulative sum, this function calculates the
1432  !! rate per unit time defined by the rate_unit argument (as default,
1433  !! per commondata::minute). This requires working on the appropriate
1434  !! difference values for each time interval.
1435  !! Note: This function version accepts *real* array (::srp).
1436  pure function cum2rate_r( inarray, rate_unit ) result(rate_out)
1437  !> Input Output_Arrays array (real ::srp type)
1438  real(srp), dimension(:), intent(in) :: inarray
1439  !> Optional rate parameter, the rate is calculated per this number of
1440  !! seconds, the default is per minute
1441  integer, intent(in), optional :: rate_unit
1442 
1443  real(srp), dimension(:), allocatable :: rate_out
1444 
1445  integer :: rate_unit_loc, oarray_max, i, j
1446 
1447  real(srp) :: sum_rate
1448 
1449  if (present(rate_unit)) then
1450  rate_unit_loc = rate_unit
1451  else
1452  rate_unit_loc = minute
1453  end if
1454 
1455  oarray_max = size(inarray) ! Output_Arrays
1456  allocate( rate_out(oarray_max/rate_unit_loc) )
1457 
1458  j = 0
1459  sum_rate = inarray(1)
1460 
1461  do i=2, oarray_max
1462  if( inarray(i) > inarray(i-1) ) then
1463  sum_rate = sum_rate + inarray(i) - inarray(i-1)
1464  end if
1465  if (mod(i,rate_unit_loc)==0) then
1466  j = j + 1
1467  rate_out(j) = sum_rate / real(rate_unit_loc,srp)
1468  sum_rate = 0.0_srp
1469  end if
1470  end do
1471 
1472  end function cum2rate_r
1473 
1474  !> Calculate rate based on simple block data, only the block start and
1475  !! end values are considered.
1476  !! @note This function is used to calculate the growth rate from body
1477  !! size dynamics data.
1478  pure function blockrate ( inarray, rate_unit ) result (rate_out)
1479  !> Input Output_Arrays array (real ::srp type)
1480  real(srp), dimension(:), intent(in) :: inarray
1481  !> Optional rate parameter, the rate is calculated per this number of
1482  !! seconds, the default is per minute
1483  integer, intent(in), optional :: rate_unit
1484 
1485  real(srp), dimension(:), allocatable :: rate_out
1486 
1487  integer :: rate_unit_loc, iarray_max, i, j
1488  real(srp) :: iarray_block_start, iarray_block_end
1489 
1490  if (present(rate_unit)) then
1491  rate_unit_loc = rate_unit
1492  else
1493  rate_unit_loc = minute
1494  end if
1495 
1496  iarray_max = size(inarray)
1497  allocate( rate_out(iarray_max/rate_unit_loc) )
1498 
1499  j=1
1500  iarray_block_start = inarray(1)
1501  do i = 2, iarray_max
1502  if ( mod(i,rate_unit_loc)==0 ) then
1503  iarray_block_end = inarray(i)
1504  rate_out(j) = (iarray_block_end - iarray_block_start) / &
1505  real(rate_unit_loc,srp)
1506  iarray_block_start = inarray(i)
1507  j=j+1
1508  end if
1509  end do
1510 
1511  end function blockrate
1512 
1513  !> Calculate an array of specific growth rate rate values based on simple
1514  !! block data,
1515  !! @note This function is identical as ::blockrate() but uses specific growth
1516  !! rate calculated by ::sgr_i() function. In turn, sgr_i()
1517  !! uses ::sgr_r() as the backend.
1518  pure function sgr_blockrate (inarray, rate_unit, time_unit) result (rate_out)
1519  !> Input Output_Arrays array (real ::srp type)
1520  real(srp), dimension(:), intent(in) :: inarray
1521  !> Optional rate parameter, the rate is calculated per this number of
1522  !! seconds, the default is per minute
1523  integer, intent(in), optional :: rate_unit
1524  !> Logical flag determining if time in s should be converted to min
1525  integer, optional, intent(in) :: time_unit
1526  real(srp), dimension(:), allocatable :: rate_out
1527 
1528  integer :: rate_unit_loc, iarray_max, i, j
1529  integer :: time_unit_loc
1530  real(srp) :: iarray_block_start, iarray_block_end
1531 
1532  if (present(rate_unit)) then
1533  rate_unit_loc = rate_unit
1534  else
1535  rate_unit_loc = minute
1536  end if
1537 
1538  if (present(time_unit)) then
1539  time_unit_loc = time_unit
1540  else
1541  time_unit_loc = 0
1542  end if
1543 
1544  iarray_max = size(inarray)
1545  allocate( rate_out(iarray_max/rate_unit_loc) )
1546 
1547  j=1
1548  iarray_block_start = inarray(1)
1549  do i = 2, iarray_max
1550  if ( mod(i,rate_unit_loc)==0 ) then
1551  iarray_block_end = inarray(i)
1552  rate_out(j) = sgr_i( iarray_block_end, iarray_block_start, &
1553  rate_unit_loc, 0, time_unit_loc)
1554  iarray_block_start = inarray(i)
1555  j=j+1
1556  end if
1557  end do
1558 
1559  end function sgr_blockrate
1560 
1561  !-----------------------------------------------------------------------------
1562  !> Simple history stack function, add to the end of the stack. We need
1563  !! only to add components on top of the stack and retain
1564  !! `commondata::history_size_spatial` elements of the prior history (for a
1565  !! spatial moving object). The stack works as follows, assuming 100 and 200
1566  !! are added:
1567  !! @n [1 2 3 4 5 6 7 8 9 10];
1568  !! @n [2 3 4 5 6 7 8 9 10 100];
1569  !! @n [3 4 5 6 7 8 9 10 100 200]
1570  !! @param history_array Integer array that keeps the history.
1571  !! @param add_this we add this element to the end of the history array.
1572  !! @note This is the integer type version.
1573  pure subroutine add_to_history_i4(history_array, add_this)
1574 
1575  ! @param history_array Integer array that keeps the history.
1576  integer, dimension(:), intent(inout) :: history_array
1577  ! @param add_this we add this element to the end of the history array.
1578  integer, intent(in) :: add_this
1579 
1580  ! Local variable to keep the size of the history array
1581  integer :: history_size
1582 
1583  history_size = size(history_array)
1584 
1585  history_array = cshift(history_array,1)
1586  history_array(history_size) = add_this
1587 
1588  end subroutine add_to_history_i4
1589 
1590  !-----------------------------------------------------------------------------
1591  !> Simple history stack function, add to the end of the stack. We need
1592  !! only to add components on top of the stack and retain
1593  !! `commondata::history_size_spatial` elements of the prior history (for a
1594  !! spatial moving object).
1595  !! @param history_array Integer array that keeps the history.
1596  !! @param add_this we add this element to the end of the history array.
1597  !! @note This is the real type version.
1598  pure subroutine add_to_history_r(history_array, add_this)
1599 
1600  ! @param history_array Integer array that keeps the history.
1601  real(srp), dimension(:), intent(inout) :: history_array
1602  ! @param add_this we add this element to the end of the history array.
1603  real(srp), intent(in) :: add_this
1604 
1605  ! Local variable to keep the size of the history array.
1606  integer :: history_size
1607 
1608  history_size = size(history_array)
1609 
1610  history_array = cshift(history_array,1)
1611  history_array(history_size) = add_this
1612 
1613  end subroutine add_to_history_r
1614 
1615  !-----------------------------------------------------------------------------
1616  !> Simple history stack function, add to the end of the stack. We need
1617  !! only to add components on top of the stack and retain
1618  !! `commondata::history_size_spatial` elements of the prior history.
1619  !! @param history_array Integer array that keeps the history.
1620  !! @param add_this we add this element to the end of the history array.
1621  !! @note This is the character string type version
1622  pure subroutine add_to_history_char(history_array, add_this)
1623 
1624  ! @param history_array Integer array that keeps the history.
1625  character(*), dimension(:), intent(inout) :: history_array
1626  ! @param add_this we add this element to the end of the history array.
1627  character(*), intent(in) :: add_this
1628 
1629  ! Local variable to keep the size of the history array.
1630  integer :: history_size
1631 
1632  history_size = size(history_array)
1633 
1634  history_array = cshift(history_array,1)
1635  history_array(history_size) = add_this
1636 
1637  end subroutine add_to_history_char
1638 
1639  !> Return the last value in the array. Note that this function is especially
1640  !! useful for long nested objects, to avoid complex `set`.
1641  pure function last_in_history_r ( history_array, offset ) result (last_value)
1642  real(srp), dimension(:), intent(in) :: history_array
1643  integer, optional, intent(in) :: offset
1644  real(srp) :: last_value
1645 
1646  integer :: offset_loc
1647 
1648  if (present(offset)) then
1649  offset_loc = offset
1650  ! offset must not go beyond the actual array limits
1651  if (offset_loc > size(history_array)-1) then
1652  offset_loc = size(history_array)-1
1653  else if (offset_loc < 0) then
1654  offset_loc = 0
1655  end if
1656  else
1657  offset_loc = 0
1658  end if
1659 
1660  last_value = history_array( size(history_array)-offset_loc )
1661 
1662  end function last_in_history_r
1663 
1664  !> Return the last value in the array. Note that this function is especially
1665  !! useful for long nested objects, to avoid complex `set`.
1666  pure function last_in_history_i4 (history_array, offset) result (last_value)
1667  integer, dimension(:), intent(in) :: history_array
1668  integer, optional, intent(in) :: offset
1669  integer :: last_value
1670 
1671  integer :: offset_loc
1672 
1673  if (present(offset)) then
1674  offset_loc = offset
1675  ! offset must not go beyond the actual array limits
1676  if (offset_loc > size(history_array)-1) then
1677  offset_loc = size(history_array)-1
1678  else if (offset_loc < 0) then
1679  offset_loc = 0
1680  end if
1681  else
1682  offset_loc = 0
1683  end if
1684 
1685  last_value = history_array( size(history_array)-offset_loc )
1686 
1687  end function last_in_history_i4
1688 
1689  !> Return the last value in the array. Note that this function is especially
1690  !! useful for long nested objects, to avoid complex `set`.
1691  pure function last_in_history_char (history_array, offset) result(last_value)
1692  character(*), dimension(:), intent(in) :: history_array
1693  integer, optional, intent(in) :: offset
1694  character(len=:), allocatable :: last_value
1695 
1696  integer :: offset_loc
1697 
1698  if (present(offset)) then
1699  offset_loc = offset
1700  ! offset must not go beyond the actual array limits
1701  if (offset_loc > size(history_array)-1) then
1702  offset_loc = size(history_array)-1
1703  else if (offset_loc < 0) then
1704  offset_loc = 0
1705  end if
1706  else
1707  offset_loc = 0
1708  end if
1709 
1710  last_value = history_array( size(history_array)-offset_loc )
1711 
1712  end function last_in_history_char
1713 
1714  !-----------------------------------------------------------------------------
1715  !> Force a value within the range set by the vmin and vmax dummy parameter
1716  !! values. If the value is within the range, it does not change, if it
1717  !! falls outside, the output force value is obtained as
1718  !! min( max( value, FORCE_MIN ), FORCE_MAX )
1719  !! @param[in] value_in Input value for forcing transformation.
1720  !! @param[in] vmin minimum value of the force-to range (lower limit), if
1721  !! not present, a lower limit of 0.0 is used.
1722  !! @param[in] vmax maximum value of the force-to range (upper limit)
1723  !! @returns an input value forced to the range.
1724  !! @note Note that this is the **real** precision version of the
1725  !! generic `within` function.
1726  elemental function within_r(value_in, vmin, vmax) result (value_out)
1727  real(srp), intent(in) :: value_in
1728  real(srp), optional, intent(in) :: vmin
1729  real(srp), intent(in) :: vmax
1730  real(srp) :: value_out
1731 
1732  ! Local copies of optionals.
1733  real(srp) :: vmin_here
1734 
1735  ! Check optional minimum value, if absent, set to a default value 0.0.
1736  if (present(vmin)) then
1737  vmin_here = vmin
1738  else
1739  vmin_here = 0.0_srp
1740  end if
1741 
1742  value_out = min( max( value_in, vmin_here ), vmax )
1743 
1744  end function within_r
1745 
1746  !-----------------------------------------------------------------------------
1747  !> Force a value within the range set by the vmin and vmax dummy parameter
1748  !! values. If the value is within the range, it does not change, if it
1749  !! falls outside, the output force value is obtained as
1750  !! min( max( value, FORCE_MIN ), FORCE_MAX )
1751  !! @param[in] value_in Input value for forcing transformation.
1752  !! @param[in] vmin minimum value of the force-to range (lower limit), if
1753  !! not present, a lower limit of 0.0 is used.
1754  !! @param[in] vmax maximum value of the force-to range (upper limit)
1755  !! @returns an input value forced to the range.
1756  !! @note Note that this is the **integer** version of the generic
1757  !! `within` function.
1758  elemental function within_i(value_in, vmin, vmax) result (value_out)
1759  integer, intent(in) :: value_in
1760  integer, optional, intent(in) :: vmin
1761  integer, intent(in) :: vmax
1762  integer :: value_out
1763 
1764  ! Local copies of optionals.
1765  integer :: vmin_here
1766 
1767  ! Check optional minimum value, if absent, set to a default value 0.0.
1768  if (present(vmin)) then
1769  vmin_here = vmin
1770  else
1771  vmin_here = 0.0_srp
1772  end if
1773 
1774  ! Calculate the force-limit output value.
1775  value_out = min( max( value_in, vmin_here ), vmax )
1776 
1777  end function within_i
1778 
1779  !> Function calculating standard oxygen consumption in rainbow trout
1780  !! based on the empirical curve defined by interpolation on the grid
1781  !! determined by two global arrays:
1782  !!
1783  !! - commondata::global_oxygen_grid_x_temp -- `smr_oxygen_temp`
1784  !! - commondata::global_oxygen_grid_y_o2std -- `smr_oxygen_o2`
1785  !! .
1786  !!
1787  !! For rainbow trout, the data are based on Evans, D.O. (1990)
1788  !! Metabolic Thermal compensation by rainbow trout: effects on
1789  !! standard metabolic rate and potential usable power.
1790  !! Trans. Am. Fish. Soc. 119, 585–600, Fig. 9:
1791  !!
1792  !! Variable |1 |2 |3 |4 |5 |6 |7
1793  !! -----------------|-----|-----|-----|-----|-----|------|-------
1794  !! Temperature ºC |0.0 |0.5 |5.1 |15.2 |20.2 |25.2 |26.2
1795  !! SMR |38.0 |38.0 |40.6 |67.8 |94.6 |136.7 |145.9
1796  !!
1797  !! Interpolation data from Evan 1990:
1798  !!@verbatim
1799  !! [0.0 0.5 5.1 15.2 20.2 25.2 26.2]
1800  !! [38.0 38.0 40.6 67.8 94.6 136.7 145.9]
1801  !!@endverbatim
1802  !!
1803  !! @note Note that this function uses `CSPLINE` cubic spline interpolation
1804  !! function from `BASE_UTILS`.
1805  pure function oxygen_rate_std_spline(temperature) result (o2rate)
1806  use base_utils, only : cspline
1807  !> Temperature Celsius. Note that temperature is an array
1808  real(srp), dimension(:), intent(in) :: temperature
1809  !> Standard Oxygen consumption @f$ mg \; O_2 \; kg ^{-1} \; h^{-1} @f$
1810  !! array for the respective input temperatures
1811  real(srp), allocatable, dimension(:) :: o2rate
1812 
1813  real(srp), dimension( size(temperature) ) :: outval
1814 
1816  temperature, outval )
1817 
1818  o2rate = outval
1819 
1820  end function oxygen_rate_std_spline
1821 
1822  !> Function calculating standard oxygen consumption in rainbow trout
1823  !! based on the empirical curve defined by interpolation on the grid
1824  !! determined by two global arrays:
1825  !!
1826  !! - commondata::global_oxygen_grid_x_temp -- `smr_oxygen_temp`
1827  !! - commondata::global_oxygen_grid_y_o2std -- `smr_oxygen_o2`
1828  !! .
1829  !!
1830  !! For rainbow trout, the data are based on Evans, D.O. (1990)
1831  !! Metabolic Thermal compensation by rainbow trout: effects on
1832  !! standard metabolic rate and potential usable power.
1833  !! Trans. Am. Fish. Soc. 119, 585–600, Fig. 9:
1834  !!
1835  !! Variable |1 |2 |3 |4 |5 |6 |7
1836  !! -----------------|-----|-----|-----|-----|-----|------|-------
1837  !! Temperature ºC |0.0 |0.5 |5.1 |15.2 |20.2 |25.2 |26.2
1838  !! SMR |38.0 |38.0 |40.6 |67.8 |94.6 |136.7 |145.9
1839  !!
1840  !! Interpolation data from Evan 1990:
1841  !!@verbatim
1842  !! [0.0 0.5 5.1 15.2 20.2 25.2 26.2]
1843  !! [38.0 38.0 40.6 67.8 94.6 136.7 145.9]
1844  !!@endverbatim
1845  !!
1846  !! @note Note that this function uses `CSPLINE_SCALAR` cubic splines
1847  !! interpolation function from `BASE_UTILS`.
1848  pure function oxygen_rate_std_ddpi(temperature) result (o2rate)
1849  use base_utils, only : cspline_scalar
1850  !> Temperature Celsius
1851  real(srp), intent(in) :: temperature
1852  !> Standard Oxygen consumption @f$ mg \; O_2 \; kg ^{-1} \; h^{-1} @f$
1853  real(srp) :: o2rate
1854 
1855  o2rate = cspline_scalar( global_oxygen_grid_x_temp, &
1856  global_oxygen_grid_y_o2std, temperature )
1857 
1858  end function oxygen_rate_std_ddpi
1859 
1860  !> Function to calculate the FCR, feed conversion rate
1861  !! @f$ FCR = \frac{F} {\Delta M} @f$
1862  elemental function fcr(input_f, output_m) result (fcr_out)
1863  !> Feed mass
1864  real(srp), intent(in) :: input_f
1865  !> Body mass increment
1866  real(srp), intent(in) :: output_m
1867  !> Feed conversion rate (FCR)
1868  real(srp) :: fcr_out
1869 
1870  fcr_out = input_f / output_m
1871 
1872  end function fcr
1873 
1874  !> Convert kg to g
1875  elemental function kg2g(kg) result (g)
1876  real(srp), intent(in) :: kg
1877  real(srp) :: g
1878 
1879  g = 1000.0_srp * kg
1880 
1881  end function kg2g
1882 
1883  !> Convert g to kg
1884  elemental function g2kg(g) result (kg)
1885  real(srp), intent(in) :: g
1886  real(srp) :: kg
1887 
1888  kg = g / 1000.0_srp
1889 
1890  end function g2kg
1891 
1892  !> Convert kg to mg
1893  elemental function kg2mg(kg) result (mg)
1894  real(srp), intent(in) :: kg
1895  real(srp) :: mg
1896 
1897  mg = 1000000.0_srp * kg
1898 
1899  end function kg2mg
1900 
1901  !> Convert mg to kg
1902  elemental function mg2kg(mg) result (kg)
1903  real(srp), intent(in) :: mg
1904  real(srp) :: kg
1905 
1906  kg = mg / 1000000.0_srp
1907 
1908  end function mg2kg
1909 
1910  !> Convert mg to g
1911  elemental function mg2g(mg) result (g)
1912  real(srp), intent(in) :: mg
1913  real(srp) :: g
1914 
1915  g = mg / 1000.0_srp
1916 
1917  end function mg2g
1918 
1919  !> Convert g to mg
1920  elemental function g2mg(g) result (mg)
1921  real(srp), intent(in) :: g
1922  real(srp) :: mg
1923 
1924  mg = g * 1000.0_srp
1925 
1926  end function g2mg
1927 
1928  !> Convert micro (\mu) unit to unit, e.g. from micro mol to mol
1929  !! Metric prefix micro = 10^-6: unit x 0.000 001
1930  elemental function from_micro(mu_unit) result (unit)
1931  real(srp), intent(in) :: mu_unit
1932  real(srp) :: unit
1933 
1934  unit = mu_unit * 0.000001_srp
1935 
1936  end function from_micro
1937 
1938  !> Calculates specific growth rate following Houde & Schekter (1981):
1939  !! @f[ G = ( e^{ \frac{ ln W_2 - ln W_1 }{ t_2 - t_1 } } -1 ) 100 @f]
1940  !! @note Note that this function is the main computation backend that
1941  !! is called in ::sgr_i() and ::sgr_blockrate().
1942  elemental function sgr_r( weight1, weight2, time1, time2 ) result (sgr_out)
1943  !> Body weight at start of the time interval
1944  real(srp), intent(in) :: weight1
1945  !> Body weight at the end of the time interval
1946  real(srp), intent(in) :: weight2
1947  !> Time at the start of the time interval
1948  real(srp), intent(in) :: time1
1949  !> Time at the end of the time interval
1950  real(srp), intent(in) :: time2
1951  real(srp) :: sgr_out
1952 
1953  sgr_out = ( exp( ( log(weight2) - log(weight1) ) / &
1954  ( time2 - time1) ) - 1 ) * 100.0_srp
1955 
1956  end function sgr_r
1957 
1958  !> Calculates specific growth rate following Houde & Schekter (1981):
1959  !! @f[ G = ( e^{ \frac{ ln W_2 - ln W_1 }{ t_2 - t_1 } } -1 ) 100 @f]
1960  !! @note Note that the `is_convert_min` (optional) parameter is useful to
1961  !! convert the time step data (s) to min.
1962  !! @note Note that this function uses ::sgr_r() as the backend.
1963  elemental function sgr_i(weight1, weight2, time1, time2, time_unit) &
1964  result(sgr_out)
1965  !> Body weight at start of the time interval
1966  real(srp), intent(in) :: weight1
1967  !> Body weight at the end of the time interval
1968  real(srp), intent(in) :: weight2
1969  !> Time at the start of the time interval
1970  integer, intent(in) :: time1
1971  !> Time at the end of the time interval
1972  integer, intent(in) :: time2
1973  !> If the two time parameters `time1` and `time2` are provided in raw
1974  !! model time units (s), this parameter allows to rescale these raw units
1975  !! to commondata::minute or commondata::hour.
1976  !! @note If time units are provided in any other measures than raw time
1977  !! steps (s), this `time_unit`parameter should be set to default
1978  !! value or out of normal range (e.g. -1). Then the rate unit will
1979  !! use the same time scale as `time1` and `time2`. Notably, this
1980  !! should ber done if the times are given in commondata::hour.
1981  integer, optional, intent(in) :: time_unit
1982  real(srp) :: sgr_out
1983 
1984  integer :: time_unit_loc
1985 
1986  if (present(time_unit)) then
1987  time_unit_loc = time_unit
1988  else
1989  time_unit_loc = 0
1990  end if
1991 
1992  if ( time_unit_loc > 1 ) then
1993  sgr_out=sgr_r(weight1, weight2, real(time1,srp)/real(time_unit_loc,srp),&
1994  real(time2,srp)/real(time_unit_loc,SRP) )
1995  else
1996  sgr_out=sgr_r(weight1, weight2, real(time1,srp), real(time2,srp) )
1997  end if
1998 
1999  end function sgr_i
2000 
2001  !> Define the default stomach emptying data structure
2002  !! ::stomach_emptying_pattern in case the input data
2003  !! from a CSV file are absent.
2004  pure function stomach_emptying_def_default() result (pattern_def)
2005  type(stomach_emptying_pattern) :: pattern_def
2006 
2007  real(srp), parameter, dimension(*) :: &
2008  mass_grid = [ 50.0_srp, 100.0_srp, 300.0_srp, 500.0_srp ]
2009 
2010  real(srp), parameter, dimension(*) :: &
2011  temp_grid = [ 5.0_srp, 10.0_srp, 15.0_srp, 20.0_srp, 22.0_srp ]
2012 
2013  integer, parameter :: levels_mass = size(mass_grid), &
2014  levels_temp = size(temp_grid)
2015 
2016  integer, parameter, dimension(LEVELS_MASS, LEVELS_TEMP) :: &
2017  emptying_grid_hour = reshape ( &
2018  ! ..........................................................................
2019  ! 50 100 300 500
2020  [ & ! 1 2 3 4 LEVELS_MASS
2021  !----------------------------
2022  24, 35, 75, 100, & ! 1 5 LEVELS_TEMP
2023  15, 20, 50, 75, & ! 2 10
2024  9, 15, 40, 50, & ! 3 15
2025  5, 10, 30, 35, & ! 4 20
2026  4, 8, 29, 33 ], & ! 5 22
2027  !-----------------------------
2028  [levels_mass, levels_temp], [0], [1,2] )
2029  ! Note: additional array shape | pad | order
2030  ! reshape params: | with | by rows
2031  ! ...........................................................................
2032 
2033  pattern_def%grid_mass = mass_grid
2034  pattern_def%grid_temperature = temp_grid
2035  pattern_def%emptying_time = emptying_grid_hour * hour
2036 
2037  end function stomach_emptying_def_default
2038 
2039  !> Show the stomach emptying data structure.
2040  impure function stomach_emptying_txt(is_raw_s) result (message_show)
2041  character(len=:), allocatable :: message_show
2042  logical, optional, intent(in) :: is_raw_s
2043 
2044  integer :: i,j
2045  logical is_raw_loc
2046  character(len=50) :: msg
2047  character(len=*), parameter :: fmt_real="f9.2", colsep=" |", &
2048  fmt_int="I9"
2049  integer :: fmt_len = 9 ! coincides with total length of FMT_REAL
2050  character(len=:), allocatable :: fmtstr, fmtint
2051 
2052  if (present(is_raw_s)) then
2053  is_raw_loc = is_raw_s
2054  else
2055  is_raw_loc = .false.
2056  end if
2057 
2058  fmtstr = "(" // fmt_real // ")"
2059  fmtint = "(" // fmt_int // ")"
2060 
2061  ! Title text
2062  message_show = "Stomach emptying time (h) parameter matrix:" // crlf
2063 
2064  ! First header line
2065  message_show = message_show // " " // colsep // &
2066  " BODY MASS:" // crlf
2067 
2068  ! Write the first row of body mass grid array
2069  message_show = message_show // " TEMP" // colsep
2070  do i=1, size( global_stomach_emptying_pattern%grid_mass )
2071  write(msg,fmtstr) global_stomach_emptying_pattern%grid_mass(i)
2072  message_show = message_show // trim(msg)
2073  end do
2074  ! Write horizontal separator line
2075  message_show = message_show // crlf // repeat("-",fmt_len) // "--+" // &
2076  repeat("-", size(global_stomach_emptying_pattern%grid_mass)*fmt_len)//crlf
2077 
2078  ! Write the main body of the data
2079  do j=1, size(global_stomach_emptying_pattern%grid_temperature)
2080  ! Each line starts with the temperature grid array value
2081  write(msg,fmtstr) global_stomach_emptying_pattern%grid_temperature(j)
2082  message_show = message_show // trim(msg) // colsep
2083  ! Normal data follow
2084  do i=1, size( global_stomach_emptying_pattern%grid_mass )
2085  if (is_raw_loc) then
2086  write(msg,fmtint) global_stomach_emptying_pattern%emptying_time(i,j)
2087  else
2088  write(msg,fmtstr) &
2089  real(global_stomach_emptying_pattern%emptying_time(i,j),srp) / hour
2090  end if
2091  message_show = message_show // trim(msg)
2092  end do
2093  message_show = message_show // crlf ! end of row
2094  end do
2095 
2096  end function stomach_emptying_txt
2097 
2098  !> Produce a formatted text string containing the adjusted stomach
2099  !! transport arrays commondata::global_transport_pattern_r and
2100  !! commondata::global_transport_pattern_t
2101  impure function stomach_transport_txt(adjust, format, show, prefix) &
2102  result(txt_out)
2103  use base_utils, only: tostr
2104  use base_strings, only: uppercase
2105  character(len=:), allocatable :: txt_out
2106  !> Optional adjustment factor for the array
2107  !! commondata::global_transport_pattern_t
2108  real(srp), optional, intent(in) :: adjust
2109  !> Optional format for the output arrays: "R" for R style, "python" for
2110  !! the Python style, default is the python style.
2111  character(len=*), optional, intent(in) :: format
2112  !> Default selection of what array to show,
2113  !! `transport_pattern_t` or `transport_pattern_r`.
2114  character(len=*), optional, intent(in) :: show
2115  !> Optional prefix string that comes before the array on left.
2116  character(len=*), optional, intent(in) :: prefix
2117 
2118  enum, bind(C)
2119  enumerator :: show_p, show_t, show_pt
2120  end enum
2121 
2122  real(srp), dimension( size(Global_Transport_Pattern_R) ) :: stpr
2123  integer, dimension( size(Global_Transport_Pattern_T) ) :: stpt
2124  real(srp) :: adjust_loc
2125  character(len=*), parameter :: fmtstr = "(f6.2)"
2126  integer :: i, what
2127  character(len=:), allocatable :: sep_l, sep_r, & ! array delimiters
2128  prefl
2129 
2130  if (present(adjust)) then
2131  adjust_loc = adjust
2132  else
2133  adjust_loc = 1.0_srp
2134  end if
2135 
2136  if (present(format)) then
2137  if (uppercase(format)=="R") then
2138  sep_l="c("
2139  sep_r=")"
2140  elseif (uppercase(format)=="PYTHON") then
2141  sep_l="["
2142  sep_r="]"
2143  end if
2144  else ! default is Python style
2145  sep_l="["
2146  sep_r="]"
2147  end if
2148 
2149  if (present(show)) then
2150  ! transport_pattern_t
2151  if (uppercase(show)=="TRANSPORT_PATTERN_T") then
2152  what=show_t
2153  ! transport_pattern_r
2154  elseif (uppercase(show)=="TRANSPORT_PATTERN_R") then
2155  what=show_p
2156  ! Both transport_pattern_t and transport_pattern_r
2157  elseif (uppercase(show)=="BOTH") then
2158  what = show_pt
2159  end if
2160  ! Default is both transport_pattern_t and transport_pattern_r
2161  else
2162  what = show_pt
2163  end if
2164 
2165  if (present(prefix)) then
2166  prefl = prefix
2167  else
2168  prefl = ""
2169  end if
2170 
2172  stpt = global_transport_pattern_t * adjust_loc
2173 
2174  if (what == show_pt .or. what==show_t) then
2175 
2176  txt_out = prefl // "transport_pattern_t = " // sep_l ! [
2177 
2178  do i=1, size(stpt)-1
2179  txt_out = txt_out // tostr(stpt(i)) // ", "
2180  end do
2181  txt_out = txt_out // tostr(stpt(size(stpt)))
2182 
2183  txt_out = txt_out // sep_r
2184 
2185  if (what==show_pt) txt_out = txt_out // crlf
2186 
2187  else
2188  txt_out = ""
2189  end if
2190 
2191  if (what==show_pt .or. what==show_p) then
2192 
2193  txt_out = txt_out // prefl // "transport_pattern_r = " // sep_l
2194 
2195  do i=1, size(stpr)-1
2196  txt_out = txt_out // tostr(stpr(i),fmtstr) // ", "
2197  end do
2198  txt_out = txt_out // tostr(stpr(size(stpr)),fmtstr)
2199 
2200  txt_out = txt_out // sep_r
2201 
2202  end if
2203 
2204  end function stomach_transport_txt
2205 
2206 end module commondata
Simple history stack function, add to the end of the stack. We need only to add components on top (en...
Definition: m_common.f90:785
Public interface for the cum2rate function that calculates rate from raw cumulative array....
Definition: m_common.f90:772
Return the last element of a history array with optional offset.
Definition: m_common.f90:792
Function calculating standard oxygen consumption in rainbow trout based on the curve from Fig....
Definition: m_common.f90:805
Function for calculating specific growth rate.
Definition: m_common.f90:811
Force a value within the range set by the vmin and vmax dummy parameter values.
Definition: m_common.f90:819
@ show_p
Definition: m_common.f90:2119
@ show_t
Definition: m_common.f90:2119
@ show_pt
Definition: m_common.f90:2119
This module defines global parameters and general-level computational utilities.
Definition: m_common.f90:13
elemental real(srp) function mg2kg(mg)
Convert mg to kg.
Definition: m_common.f90:1903
character(len=1), parameter, public prog_bar_char
The symbol used for writing the progress bar in text mode.
Definition: m_common.f90:77
character(len=file_name_len), public global_food_pattern_file
Global variable that keeps the file name for the food provisioning pattern see environ::food_provisio...
Definition: m_common.f90:641
pure logical function is_stress()
Calculate if stress factor grid is enabled. Stress is disabled if any of the stress parameters are no...
Definition: m_common.f90:1196
logical, public global_stress_intervention_is_minutes
Global logical flag setting the time unit for the stress intervention parameter array defined by glob...
Definition: m_common.f90:392
character(len= *), parameter, public output_dest_dir
Default output directory: all plot and data file will be saved into this directory by default (local ...
Definition: m_common.f90:99
integer, dimension(:), allocatable, public global_transport_pattern_t
Food transport in stomach: Time grid array for interpolation. Configuration file example:
Definition: m_common.f90:325
@, public unscaled
Definition: m_common.f90:177
real(srp), public global_food_input_rate
Standard rate of food item input, per minute. Configuration file example:
Definition: m_common.f90:618
pure real(srp) function, dimension(:), allocatable blockrate(inarray, rate_unit)
Calculate rate based on simple block data, only the block start and end values are considered.
Definition: m_common.f90:1479
elemental real(srp) function g2kg(g)
Convert g to kg.
Definition: m_common.f90:1885
elemental real(srp) function midgut_temp_factor(temperature)
Calculate temperature adjustment factor for the Michaelis-Menten absorption process in the midgut....
Definition: m_common.f90:1255
character(len= *), parameter, public svn_version_global
Definition: m_common.f90:31
real(srp), public global_food_item_mass
Standard mass of one food item. Configuration file example:
Definition: m_common.f90:593
integer, public global_digestion_delay_min
Delay of digestion, min. It is the delay of absorption after a food item was transmitted to the mid-g...
Definition: m_common.f90:275
real(srp), public global_feeding_activity_factor
Locomotor activity multiplier factor for the activity during the feeding time, i.e....
Definition: m_common.f90:690
real(srp), public global_water_uptake_a
Parameters of the logistic water uptake function that defines the temporary pattern of water uptake c...
Definition: m_common.f90:303
impure character(len=:) function, allocatable stomach_transport_txt(adjust, format, show, prefix)
Produce a formatted text string containing the adjusted stomach transport arrays commondata::global_t...
Definition: m_common.f90:2103
real(srp), public global_baseline_activity_day
Global baseline locomotor activity pattern (swimming speed) for each time step consists of two compon...
Definition: m_common.f90:685
integer, dimension(2), public global_interval_food_param
Parameters of the food provisioning pattern:
Definition: m_common.f90:629
elemental real(srp) function o2vol(mass_g)
Calculate the volume of 1 dm3 (l) of oxygen.
Definition: m_common.f90:1351
logical, public global_food_pattern_file_is_propagate
Global logical flag that defines if the feed scheduling pattern defined by the global_interval_food_p...
Definition: m_common.f90:649
elemental real(srp) function, private sgr_i(weight1, weight2, time1, time2, time_unit)
Calculates specific growth rate following Houde & Schekter (1981):
Definition: m_common.f90:1965
character(len= *), parameter, public def_stats_output_file
Default value of the general simulation statistics file defined in commondata::global_stats_output_fi...
Definition: m_common.f90:722
real(srp), public global_midgut_mass
Fish mid-gut mass capacity, max. filling capacity. Configuration file example:
Definition: m_common.f90:240
integer, public global_day_starts_hour
Default hour at which the "daytime" is normally started. This means, in particular,...
Definition: m_common.f90:204
elemental integer function, private within_i(value_in, vmin, vmax)
Force a value within the range set by the vmin and vmax dummy parameter values. If the value is withi...
Definition: m_common.f90:1759
real(srp), public global_appetite_fish_night
Fish appetite at night. This value is normally low because the fish do not feed at night....
Definition: m_common.f90:459
real(srp), public global_sda_absorp_rate_max
Defines the specific dynamic action (SDA) that depends on the absorption rate. There are two paramete...
Definition: m_common.f90:554
character(len= *), parameter, public output_arrays_csv_file
Default file base-name for saving model output arrays.
Definition: m_common.f90:89
real(srp), dimension(:), allocatable, public global_oxygen_grid_y_o2std
This parameter defines the Y axis of the grid: oxygen consumption. See commondata::oxygen_rate_std() ...
Definition: m_common.f90:586
logical, public global_stomach_midgut_mass_is_automatic
Logical flag that specifies that the stomach and midgut mass (capacity) are automatically recalculate...
Definition: m_common.f90:249
integer, public global_transport_dim
Food transport in stomach: Dimensionality of the transport pattern in the stomach.
Definition: m_common.f90:319
integer, public global_run_model_hours
Default duration of the model run in hours. Note that parameter file uses hours as unit....
Definition: m_common.f90:196
pure real(srp) function, private oxygen_rate_std_ddpi(temperature)
Function calculating standard oxygen consumption in rainbow trout based on the empirical curve define...
Definition: m_common.f90:1849
real(srp), public global_body_mass
Fish body mass, g, at the start of the simulation Configuration file example:
Definition: m_common.f90:226
type(stomach_emptying_pattern), public global_stomach_emptying_pattern
Default stomach emptying pattern.
Definition: m_common.f90:451
elemental real(srp) function michaelis_menten(x, r_max, x50)
Michaelis-Meneten function:
Definition: m_common.f90:1236
logical, parameter, public is_debug
Logical flag that sets the DEBUG mode. When the program is running in the DEBUG mode,...
Definition: m_common.f90:40
real(srp), public global_mid_gut_mm_k_m
Michaelis-Meneten food absorption parameter in mid-gut, (see the_fish::michaelis_menten()) relative ...
Definition: m_common.f90:499
character(len= *), parameter, public output_rate_data_csv_file
Default file base-name for saving model rate data.
Definition: m_common.f90:93
pure real(srp) function, dimension(:), allocatable cum2rate_int(inarray, rate_unit)
Given an in put array of cumulative sum, this function calculates the rate per unit time defined by t...
Definition: m_common.f90:1396
elemental real(srp) function, private within_r(value_in, vmin, vmax)
Force a value within the range set by the vmin and vmax dummy parameter values. If the value is withi...
Definition: m_common.f90:1727
real(srp), public global_ue_ze_ammonia_excretion
A parameter defining branchial and urine (ZE+UE) energy consumption as a fixed ammonia excretion rate...
Definition: m_common.f90:535
character(len=file_name_len), public global_param_file_name
Defines the name of the file that keeps global parameters. This name is normally defined by the envir...
Definition: m_common.f90:153
integer, parameter, public file_name_len
The default length of the file name.
Definition: m_common.f90:54
real(srp), public global_stress_activity_decr
Maximum value of the suppressive effect of stress on baseline activity. The actual suppression of the...
Definition: m_common.f90:383
elemental real(srp) function energy_o2(l_o2)
Calculate the enegretic equivalent (kJ) of Oxygen uptake (l O2) Source: Rønnestad B291.
Definition: m_common.f90:1326
elemental real(srp) function kg2g(kg)
Convert kg to g.
Definition: m_common.f90:1876
integer, public global_run_model_feed_offset
The offset (delay) to start feeding at the beginning of the simulation, Note that the parameter file ...
Definition: m_common.f90:211
real(srp), public global_stomach_mass
Fish stomach mass, max. filling capacity. Configuration file example:
Definition: m_common.f90:233
real(srp), public global_baseline_activity_night
Definition: m_common.f90:685
impure character(len=:) function, allocatable model_parameters_txt()
Produce a text string containing all global parameters of the model.
Definition: m_common.f90:835
logical, parameter, public is_progress_bar_gui
Global parameter defining if progress bar is printed in the GUI mode.
Definition: m_common.f90:68
integer, parameter, public minute
Definition: m_common.f90:61
pure subroutine, private add_to_history_char(history_array, add_this)
Simple history stack function, add to the end of the stack. We need only to add components on top of ...
Definition: m_common.f90:1623
integer, public global_maximum_duration_midgut_min
The maximum duration a food item can be processed in the fish mid-gut, min. If it stays in the mid-gu...
Definition: m_common.f90:284
character(len=file_name_len), public global_stats_output_file
Global variable keeping the file name for saving general simulation statistics. Data (rows) for each ...
Definition: m_common.f90:699
real(srp), public global_stress_cost_smr
Maximum value of the metabolic cost of stress in units of resting metabolic rate (SMR)....
Definition: m_common.f90:374
real(srp), public global_absorption_ratio
Maximum absorption ratio relative to the original dry food item mass. Configuration file example:
Definition: m_common.f90:257
elemental real(srp) function from_micro(mu_unit)
Convert micro (\mu) unit to unit, e.g. from micro mol to mol Metric prefix micro = 10^-6: unit x 0....
Definition: m_common.f90:1931
real(srp), dimension(*), parameter trout_oxygen_grid_x_temp
Interpolation grid defining the function calculating standard oxygen consumption in rainbow trout bas...
Definition: m_common.f90:743
character(len=file_name_len), public global_output_dest
Output directory: all plot and output files will be saved into this directory, default value is set b...
Definition: m_common.f90:134
logical, parameter, public is_stress_enable
Compile-time logical flag that defines if stress effects described by global_stress_factor_hour,...
Definition: m_common.f90:45
integer, public global_ingestion_delay_min
Delay of ingestion, min.
Definition: m_common.f90:266
real(srp), dimension(:), allocatable, public global_temp_factor_midgut_m
Temperature adjustment for absorption is controlled by the two parameters below that define the tempe...
Definition: m_common.f90:415
pure real(srp) function, dimension(:), allocatable, private sgr_blockrate(inarray, rate_unit, time_unit)
Calculate an array of specific growth rate rate values based on simple block data,...
Definition: m_common.f90:1519
integer, parameter, public hour
Global time constants, number of sec in hour and min.
Definition: m_common.f90:61
pure real(srp) function, private last_in_history_r(history_array, offset)
Return the last value in the array. Note that this function is especially useful for long nested obje...
Definition: m_common.f90:1642
@, public scale_min
Definition: m_common.f90:177
character(len= *), parameter, public param_file_name_def
Defines the name of the file that keeps global parameters. The parameter file contains the coupled pa...
Definition: m_common.f90:146
pure subroutine, private add_to_history_r(history_array, add_this)
Simple history stack function, add to the end of the stack. We need only to add components on top of ...
Definition: m_common.f90:1599
real(srp), public global_appetite_stomach_threshold
Protective appetite threshold for stomach: this is the maximum value of the the_fish::stomach::appeti...
Definition: m_common.f90:469
real(srp), dimension(:), allocatable, public global_stress_fact_suppress
The pattern of stress effect on the appetite is described by two grid arrays global_stress_factor_hou...
Definition: m_common.f90:362
character(len= *), parameter, public txt_sec
Time scale abbreviation strings for Y axes of output plots.
Definition: m_common.f90:764
elemental real(srp) function mg2g(mg)
Convert mg to g.
Definition: m_common.f90:1912
real(srp), public global_transport_baseline_fish_mass
The baseline fish mass that applies to the stomach transport pattern defined by global_transport_patt...
Definition: m_common.f90:429
character(len=file_name_len), public global_stomach_emptying_matrix_file
Global variable keeping the file name for the baseline stomach emptying matrix file that keeps the st...
Definition: m_common.f90:710
elemental real(srp) function feed_energy(feed_mass, gross_energy, is_adjust_water_uptake, water_uptake)
Calculate the energy content of a feed mass.
Definition: m_common.f90:1277
real(srp), public global_food_gross_energy
Gross energy content of the feed, MJ/kg (=kJ7g) Configuration file example:
Definition: m_common.f90:600
@, public scale_sec
Definition: m_common.f90:177
real(srp), parameter, public base_e
Natural logarithm base, e
Definition: m_common.f90:725
logical, public global_output_stats_is_long
Define if the output stats are saved using the long or short format. In the first case,...
Definition: m_common.f90:718
character(len= *), parameter, public program_title
Program title that appears on the main interface/window.
Definition: m_common.f90:34
real(srp), dimension(:), allocatable, public global_stress_factor_hour
The pattern of stress effect on the appetite is described by two grid arrays global_stress_factor_hou...
Definition: m_common.f90:352
integer, parameter, public gui_progress_bar_steps
Number of steps in the progress bar at the GUI mode.
Definition: m_common.f90:74
elemental real(srp) function kg2mg(kg)
Convert kg to mg.
Definition: m_common.f90:1894
real(srp), dimension(:), allocatable, public global_temp_factor_midgut_t
Temperature adjustment for absorption is controlled by the two parameters below that define the tempe...
Definition: m_common.f90:342
@, public scale_hour
Definition: m_common.f90:177
integer, parameter, public label_len
The default length of labels and similar text strings.
Definition: m_common.f90:51
logical, public global_is_ue_ze_fixed_rate
Logical flag to choose how branchial and urinal (ZE+UE) energy loss is calculated.
Definition: m_common.f90:544
real(srp), parameter, public debug_warn_level
Warning level in the debugging mode.
Definition: m_common.f90:48
integer, parameter, public minimum_filename_chars
The minimum length of a data file. Note that cancelled file name dialog returns a six-char random seq...
Definition: m_common.f90:58
real(srp), public global_temperature
Ambient temperature Configuration file example:
Definition: m_common.f90:611
real(srp), parameter, public min_mass_toler
Minimum mass of food item that counts as non-zero, a tolerance limit.
Definition: m_common.f90:604
real(srp), public global_sda_factor_max
Defines the specific dynamic action (SDA) that depends on the absorption rate. There are two paramete...
Definition: m_common.f90:566
integer, parameter, public prog_bar_len
The length of the progress bar in text mode.
Definition: m_common.f90:80
elemental real(srp) function fcr(input_f, output_m)
Function to calculate the FCR, feed conversion rate .
Definition: m_common.f90:1863
integer, public global_rate_interval
Default rate discretization interval in minutes, this interval is used to plot the ingestion rate.
Definition: m_common.f90:188
elemental real(srp) function g2mg(g)
Convert g to mg.
Definition: m_common.f90:1921
real(srp), public global_mid_gut_mm_r_max
Michaelis-Meneten food absorption parameter in mid-gut, (see the_fish::michaelis_menten()) relative ...
Definition: m_common.f90:491
integer(long), dimension(:), allocatable, public global_stress_intervention_time
The time of stress intervention is set by this array. It is the start points of stressful effects (s)...
Definition: m_common.f90:403
pure real(srp) function, dimension(:), allocatable, private oxygen_rate_std_spline(temperature)
Function calculating standard oxygen consumption in rainbow trout based on the empirical curve define...
Definition: m_common.f90:1806
real(srp), public global_appetite_activity_factor
Activity appetite factor determining how fish locomotor activity increases with increasing appetite....
Definition: m_common.f90:519
pure subroutine, private add_to_history_i4(history_array, add_this)
Simple history stack function, add to the end of the stack. We need only to add components on top of ...
Definition: m_common.f90:1574
real(srp), public global_energy_appetite_rate
The steepness parameter of the Logistic energy component of appetite.
Definition: m_common.f90:505
integer(long), public global_time_step
Global variable that defines the current time step.
Definition: m_common.f90:86
integer, parameter, public gui_progress_bar_max
Maximum value (steps) in the progress bar (minimum is 0)
Definition: m_common.f90:71
character(len= *), parameter, public txt_hour
Definition: m_common.f90:764
pure real(srp) function, dimension(:), allocatable cum2rate_r(inarray, rate_unit)
Given an in put array of cumulative sum, this function calculates the rate per unit time defined by t...
Definition: m_common.f90:1437
real(srp), public global_appetite_logist_r
Logistic function parameter R for the appetite factor. See the_fish::appetite_func().
Definition: m_common.f90:483
logical, dimension(:), allocatable, public global_interval_food_pattern
Global food provisioning pattern, logical TRUE/FALSE for each time step.
Definition: m_common.f90:632
character, parameter, public crlf
Definition: m_common.f90:83
real(srp), public global_energy_appetite_shift
The shift parameter of the Logistic energy component of appetite.
Definition: m_common.f90:511
real(srp), public global_appetite_logist_a
Logistic function parameter A for the appetite factor. See the_fish::appetite_func().
Definition: m_common.f90:476
pure integer function, private last_in_history_i4(history_array, offset)
Return the last value in the array. Note that this function is especially useful for long nested obje...
Definition: m_common.f90:1667
elemental integer(long) function rate2int(rate_min)
The function converts rate per min to interval in seconds.
Definition: m_common.f90:1363
pure character(len=:) function, allocatable, private last_in_history_char(history_array, offset)
Return the last value in the array. Note that this function is especially useful for long nested obje...
Definition: m_common.f90:1692
elemental real(srp) function, private sgr_r(weight1, weight2, time1, time2)
Calculates specific growth rate following Houde & Schekter (1981):
Definition: m_common.f90:1943
logical, public verbose
Global logical flag to suppress extra text diagnostics and messages, the "quiet mode" or the "verbose...
Definition: m_common.f90:173
pure type(stomach_emptying_pattern) function stomach_emptying_def_default()
Define the default stomach emptying data structure ::stomach_emptying_pattern in case the input data ...
Definition: m_common.f90:2005
real(srp), dimension(:), allocatable, public global_transport_pattern_r
Food transport in stomach: Proportion of food mass left in stomach. Configuration file example:
Definition: m_common.f90:331
integer, public global_hours_daytime_feeding
Default duration of the day time in hours, night duration is defined as 24 - ::global_hours_day_feedi...
Definition: m_common.f90:219
real(srp), public global_water_uptake_r
Parameters of the logistic water uptake function that defines the temporary pattern of water uptake c...
Definition: m_common.f90:315
real(srp), public global_ue_ze_factor
A parameter defining branchial and urine (ZE+UE) energy consumption as a factor to SMR,...
Definition: m_common.f90:527
impure character(len=:) function, allocatable stomach_emptying_txt(is_raw_s)
Show the stomach emptying data structure.
Definition: m_common.f90:2041
real(srp), dimension(:), allocatable, public global_oxygen_grid_x_temp
This parameter defines the X axis of the grid: temperature. See commondata::oxygen_rate_std() for det...
Definition: m_common.f90:576
character(len= *), parameter, public txt_min
Definition: m_common.f90:764
elemental real(srp) function o2mass(vol_l)
Calculate the mass of 1 dm3 (l) of Oxygen.
Definition: m_common.f90:1339
integer, parameter, public max_food_items_index
Maximum number of food items in the fish stomach, this is the maximum index for the food items array ...
Definition: m_common.f90:168
real(srp), public global_transport_baseline_temperature
The baseline temperature that applies to the stomach transport pattern defined by global_transport_pa...
Definition: m_common.f90:422
elemental real(srp) function logistic(X, A, C, R)
Logistic function is defined by the equation.
Definition: m_common.f90:1214
integer, parameter, public max_points_plot
The maximum number of data points that are plotted for the model output arrays (squeezed arrays).
Definition: m_common.f90:65
real(srp), dimension(*), parameter trout_oxygen_grid_y_o2std
Interpolation grid defining the function calculating standard oxygen consumption in rainbow trout bas...
Definition: m_common.f90:760
logical, public global_food_pattern_file_is_steps
Global logical flag that defines if the feed scheduling pattern defined by the global_food_pattern_fi...
Definition: m_common.f90:658
real(srp), public global_water_uptake
Proportion of water uptake, relative of dry mass of the food item. Configuration file example:
Definition: m_common.f90:291
A data structure that defines the stomach emptying pattern as dependent on the temperature and the fi...
Definition: m_common.f90:444