Line data Source code
1 : !> \file mo_mpr_eval.f90
2 : !> \brief \copybrief mo_mpr_eval
3 : !> \details \copydetails mo_mpr_eval
4 :
5 : !> \brief Runs MPR
6 : !> \details Runs MPR and writes to global effective parameters
7 : !> \authors Robert Schweppe
8 : !> \date Feb 2018
9 : !> \copyright Copyright 2005-\today, the mHM Developers, Luis Samaniego, Sabine Attinger: All rights reserved.
10 : !! mHM is released under the LGPLv3+ license \license_note
11 : !> \ingroup f_mpr
12 : MODULE mo_mpr_eval
13 :
14 : USE mo_kind, ONLY : i4, dp
15 :
16 : IMPLICIT NONE
17 :
18 : PRIVATE
19 :
20 : PUBLIC :: mpr_eval
21 :
22 : ! ------------------------------------------------------------------
23 :
24 : CONTAINS
25 :
26 : ! ------------------------------------------------------------------
27 :
28 : ! NAME
29 : ! mpr_eval
30 :
31 : ! PURPOSE
32 : !> \brief Runs MPR and writes to global effective parameters
33 :
34 : !> \details Runs MPR and writes to global effective parameters
35 :
36 : ! INTENT(IN), OPTIONAL
37 : !> \param[in] "real(dp), dimension(:), optional :: parameterset" a set of global parameter (gamma) to run mHM,
38 : !> DIMENSION [no. of global_Parameters]
39 :
40 : ! HISTORY
41 : !> \authors Juliane Mai, Rohini Kumar
42 :
43 : !> \date Feb 2013
44 :
45 : ! Modifications:
46 : ! R. Kumar Jun 2013 - restart_flag_states_read is passed to mhm call for the soil moisture initalisation
47 : ! R. Kumar Jun 2013 - frac_sealed_city_area is added
48 : ! R. Kumar & S. Thober Aug 2013 - code change to incorporate output timestep during writing of the netcdf file
49 : ! R. Kumar Aug 2013 - added iFlag_LAI_data_format to handle LAI options, and changed within the code made accordingly
50 : ! R. Kumar, J. Mai Sep 2013 - Splitting allocation and initialization of arrays
51 : ! R. Kumar Nov 2013 - update intent variables in documentation
52 : ! L. Samaniego Nov 2013 - relational statements == to .eq., etc.
53 : ! M. Zink Feb 2014 - added PET calculation: Hargreaves-Samani (Process 5)
54 : ! M. Zink Mar 2014 - added inflow from upstream areas
55 : ! Stephan Thober Jun 2014 - added chunk read for meteorological input
56 : ! Stephan Thober Jun 2014 - updated flag for read_restart
57 : ! M. Cuntz & J. Mai Nov 2014 - LAI input from daily, monthly or yearly files
58 : ! Matthias Zink Dec 2014 - adopted inflow gauges to ignore headwater cells
59 : ! Stephan Thober Aug 2015 - moved writing of daily discharge to mo_write_routing, included routing related variables from mRM
60 : ! David Schaefer Aug 2015 - changed to new netcdf-writing scheme
61 : ! Stephan Thober Sep 2015 - updated mrm_routing call
62 : ! O. Rakovec, R. Kumar Oct 2015 - added optional output for domain averaged TWS
63 : ! Rohini Kumar Mar 2016 - changes for handling multiple soil database options
64 : ! Stephan Thober Nov 2016 - added two options for routing
65 : ! Rohini Kuamr Dec 2016 - option to handle monthly mean gridded fields of LAI
66 : ! Stephan Thober Jan 2017 - added prescribed weights for tavg and pet
67 : ! Zink M. Demirel C. Mar 2017 - Added Jarvis soil water stress function at SM process(3)
68 : ! Robert Schweppe Dec 2017 - extracted call to mpr from inside mhm
69 :
70 60 : SUBROUTINE mpr_eval(parameterset, opti_domain_indices)
71 :
72 : use mo_common_variables, only : L0_LCover, l0_l1_remap, level0, level1, domainMeta
73 : use mo_message, only : message
74 : use mo_mpr_global_variables, only : L0_asp, L0_geoUnit, L0_gridded_LAI, &
75 : L0_slope_emp, L0_soilId, L1_HarSamCoeff, L1_PrieTayAlpha, L1_aeroResist, &
76 : L1_alpha, L1_degDayInc, L1_degDayMax, L1_degDayNoPre, L1_fAsp, L1_fRoots, &
77 : L1_fSealed, L1_jarvis_thresh_c1, L1_kBaseFlow, L1_kPerco, L1_kSlowFlow, &
78 : L1_karstLoss, L1_kfastFlow, L1_maxInter, L1_petLAIcorFactor, L1_sealedThresh, &
79 : L1_soilMoistExp, L1_soilMoistFC, L1_soilMoistSat, L1_surfResist, L1_tempThresh, &
80 : L1_unsatThresh, L1_wiltingPoint, &
81 : ! neutron count
82 : L1_No_Count, L1_bulkDens, L1_latticeWater, L1_COSMICL3
83 :
84 : use mo_multi_param_reg, only : mpr
85 : use mo_string_utils, only : num2str
86 : use mo_timer, only : timer_clear, timer_get, timer_start, timer_stop
87 :
88 : implicit none
89 :
90 : ! a set of global parameter (gamma) to run mHM, DIMENSION [no. of global_Parameters]
91 : real(dp), dimension(:), intent(in), optional :: parameterset
92 :
93 : integer(i4), dimension(:), optional, intent(in) :: opti_domain_indices
94 :
95 : ! Counters
96 : integer(i4) :: iDomain, ii, nDomains, itimer
97 :
98 : ! start and end index at level 0 for current domain
99 : integer(i4) :: s0, e0
100 :
101 : ! start and end index at level 1 for current domain
102 : integer(i4) :: s1, e1
103 :
104 :
105 : !-------------------------------------------------------------------
106 : ! NOW call MPR
107 : !-------------------------------------------------------------------
108 60 : call message(' Executing MPR ...')
109 60 : itimer = 10
110 60 : call timer_start(itimer)
111 :
112 : !----------------------------------------
113 : ! loop over domains
114 : !----------------------------------------
115 60 : if (present(opti_domain_indices)) then
116 0 : nDomains = size(opti_domain_indices)
117 : else
118 60 : nDomains = domainMeta%nDomains
119 : end if
120 267 : do ii = 1, nDomains
121 207 : if (present(opti_domain_indices)) then
122 0 : iDomain = opti_domain_indices(ii)
123 : else
124 : iDomain = ii
125 : end if
126 :
127 : ! get domain information
128 207 : s0 = level0(domainMeta%L0DataFrom(iDomain))%iStart
129 207 : e0 = level0(domainMeta%L0DataFrom(iDomain))%iEnd
130 207 : s1 = level1(iDomain)%iStart
131 207 : e1 = level1(iDomain)%iEnd
132 :
133 0 : call mpr(level0(domainMeta%L0DataFrom(iDomain))%mask, L0_geoUnit(s0 : e0), &
134 0 : L0_soilId(s0 : e0, :), L0_asp(s0 : e0), L0_gridded_LAI(s0 : e0, :), &
135 0 : L0_LCover(s0 : e0, :), L0_slope_emp(s0 : e0), &
136 207 : pack(level0(domainMeta%L0DataFrom(iDomain))%y, level0(domainMeta%L0DataFrom(iDomain))%mask), &
137 0 : level0(domainMeta%L0DataFrom(iDomain))%Id, &
138 0 : l0_l1_remap(iDomain)%upper_bound, l0_l1_remap(iDomain)%lower_bound, &
139 0 : l0_l1_remap(iDomain)%left_bound, l0_l1_remap(iDomain)%right_bound, &
140 0 : l0_l1_remap(iDomain)%n_subcells, &
141 0 : L1_fSealed(s1 : e1, :, :), &
142 0 : L1_alpha(s1 : e1, :, :), L1_degDayInc(s1 : e1, :, :), L1_degDayMax(s1 : e1, :, :), &
143 0 : L1_degDayNoPre(s1 : e1, :, :), L1_fAsp(s1 : e1, :, :), L1_HarSamCoeff(s1 : e1, :, :), &
144 0 : L1_PrieTayAlpha(s1 : e1, :, :), L1_aeroResist(s1 : e1, :, :), L1_surfResist(s1 : e1, :, :), &
145 0 : L1_fRoots(s1 : e1, :, :), L1_kFastFlow(s1 : e1, :, :), &
146 0 : L1_kSlowFlow(s1 : e1, :, :), L1_kBaseFlow(s1 : e1, :, :), L1_kPerco(s1 : e1, :, :), &
147 0 : L1_karstLoss(s1 : e1, :, :), L1_soilMoistFC(s1 : e1, :, :), L1_soilMoistSat(s1 : e1, :, :), &
148 0 : L1_soilMoistExp(s1 : e1, :, :), L1_jarvis_thresh_c1(s1 : e1, :, :), &
149 0 : L1_tempThresh(s1 : e1, :, :), L1_unsatThresh(s1 : e1, :, :), L1_sealedThresh(s1 : e1, :, :), &
150 0 : L1_wiltingPoint(s1 : e1, :, :), L1_maxInter(s1 : e1, :, :), L1_petLAIcorFactor(s1 : e1, :, :), &
151 0 : L1_No_Count(s1:e1,:,:), L1_bulkDens(s1:e1,:,:), L1_latticeWater(s1:e1,:,:), L1_COSMICL3(s1:e1,:,:), &
152 474 : parameterset )
153 :
154 : end do
155 60 : call timer_stop(itimer)
156 60 : call message(' in ', trim(num2str(timer_get(itimer), '(F9.3)')), ' seconds.')
157 60 : call timer_clear(itimer)
158 :
159 60 : END SUBROUTINE mpr_eval
160 :
161 : END MODULE mo_mpr_eval
|