LandscapeDNDC 1.37.0
Loading...
Searching...
No Matches
Public Member Functions | Private Member Functions | Private Attributes | List of all members
ldndc::EcHy Class Reference

Hydrology model EcHy. More...

#include <models/watercycle/echy/echy.h>

Public Member Functions

lerr_t solve ()
 

Private Member Functions

lerr_t EcHyIrrigation ()
 Irrigation.
 
lerr_t EcHyFlood ()
 sets hydrologic conditions during flooding events, e.g.,
 
lerr_t EcHyPercolation (size_t, size_t)
 Calculates water percolation within the soil profile.
 
lerr_t EcHyBypassFlow (double, double)
 Calculates water percolation within the soil profile.
 
lerr_t EcHyEvapotranspiration ()
 Calculates evapotranspiration within the soil profile.
 
lerr_t EcHyCapillaryRise (double const &, double const &, size_t const &, double const &, double &)
 ...
 
lerr_t EcHyBalanceCheck (double &)
 ...
 
void EcHyreset ()
 ...
 
lerr_t EcHyGroundwater ()
 ...
 
lerr_t EcHyIntegration ()
 ...
 
lerr_t EcHyStepInit ()
 ...
 
lerr_t EcHyStepExit ()
 ...
 
double EcHyGetInterceptionCapacity ()
 ...
 
lerr_t EcHyCalculateLeafWaterDistribution (lvector_t< double > &, lvector_t< double > &)
 ...
 
lerr_t EcHyPotentialEvapotranspiration ()
 Calculates potential evapotranspiration. Specific concept can be given as model option.
 
lerr_t EcHySnowIce ()
 Calls SnowDNDC for the calculation of snowpack and soil ice formation.
 
double EcHyGetMinimumWater (size_t)
 ...
 
double EcHyGetWiltingPoint (size_t)
 ...
 
double EcHyGetAvailableWaterTranspiration (size_t)
 ...
 
double EcHyGetAvailableWaterEvaporation (size_t)
 ...
 
double EcHyGetRootLimitation (size_t _sl, double)
 ...
 
double EcHySoilWaterChange (size_t)
 

Private Attributes

double kst_bottom
 
double gw_depth_static
 
double accumulated_potentialtranspiration_old
 
double accumulated_irrigation_old
 
lvector< double > trwl_sl
 
lvector< double > wlfl_sl
 
lvector< double > cr_fill_sl
 
lvector< double > gw_fill_sl
 
lvector< double > bypass_fill_sl
 
double ev_surfacewater
 
lvector< double > evsws_sl
 
lvector< double > kst_sl
 
lvector< double > wl_sl
 
lvector< double > wlfc_sl
 
lvector< double > wlwp_sl
 
double ev_leaf
 
double gw_fill_surface
 
double cr_fill_groundwater
 
double runoff
 
double throughfall_snow
 
double throughfall_water
 
double canopy_snowmelt
 
double canopy_water
 
double thornthwaite_heat_index
 
double daily_potential_evapotranspiration
 
double daily_potential_leaf_evaporation
 
double daily_potential_soil_evaporation
 
double daily_potential_transpiration
 
WaterCycleSnowDNDC::IceContentStateIn m_icecontent_in
 
cbm::string_t evapotranspiration_method
 
cbm::string_t runoff_method
 

Detailed Description

Hydrology model EcHy.

Author
David Kraus

Member Function Documentation

◆ EcHyBalanceCheck()

lerr_t ldndc::EcHy::EcHyBalanceCheck ( double & _balance)
private

...

Checks balance between all incoming and outgoing water fluxes.

570{
571 double balance( canopy_water + surface_water + canopy_snow + surface_snow);
572 for ( size_t sl = 0; sl < soillayers_in->soil_layer_cnt(); ++sl)
573 {
574 balance += wl_sl[sl] + wc_.ice_sl[sl] * sc_.h_sl[sl];
575 }
576
577 if ( _balance > 0.0)
578 {
579 balance += wlfl_sl[soillayers_in->soil_layer_cnt()]
580 + ev_leaf + ev_surfacewater + evsws_sl.sum() + trwl_sl.sum() + runoff
581 - gw_fill_surface - cr_fill_groundwater - gw_fill_sl.sum()
582 - irrigation - throughfall_water - throughfall_snow;
583
584 double const balance_delta( std::abs( _balance - balance));
585 if ( cbm::flt_greater( balance_delta, 1.0e-4))
586 {
587 KLOGWARN( " Water leakage in: ", name(),
588 " Difference: ", balance - _balance);
589 return LDNDC_ERR_FAIL;
590 }
591 }
592 else
593 {
594 _balance = balance;
595 }
596
597 return LDNDC_ERR_OK;
598}
ldndc::EcHy::gw_fill_surface
double gw_fill_surface
Definition echy.h:129
ldndc::EcHy::wl_sl
lvector< double > wl_sl
Definition echy.h:117
ldndc::EcHy::runoff
double runoff
Definition echy.h:135
ldndc::EcHy::ev_surfacewater
double ev_surfacewater
Definition echy.h:107
ldndc::EcHy::throughfall_water
double throughfall_water
Definition echy.h:141
ldndc::EcHy::cr_fill_groundwater
double cr_fill_groundwater
Definition echy.h:132
ldndc::EcHy::evsws_sl
lvector< double > evsws_sl
Definition echy.h:110
ldndc::EcHy::throughfall_snow
double throughfall_snow
Definition echy.h:138
ldndc::EcHy::trwl_sl
lvector< double > trwl_sl
Definition echy.h:85
ldndc::EcHy::canopy_water
double canopy_water
Definition echy.h:148
ldndc::EcHy::ev_leaf
double ev_leaf
Definition echy.h:126
ldndc::EcHy::gw_fill_sl
lvector< double > gw_fill_sl
Definition echy.h:101
ldndc::EcHy::wlfl_sl
lvector< double > wlfl_sl
Definition echy.h:93

References canopy_water, cr_fill_groundwater, ev_leaf, ev_surfacewater, evsws_sl, gw_fill_sl, gw_fill_surface, runoff, throughfall_snow, throughfall_water, trwl_sl, wl_sl, and wlfl_sl.

Referenced by EcHyEvapotranspiration(), and solve().

Here is the caller graph for this function:

◆ EcHyCapillaryRise()

lerr_t ldndc::EcHy::EcHyCapillaryRise ( double const & _cp,
double const & _delta_z,
size_t const & _sl,
double const & _saturation,
double & _flow )
private

...

This routine calculates the rate of capillary flow or percolation between groundwater table and root zone. The stationary flow is found by integration of dZL = K.d(MH)/(K + FLW), where Z= height above groundwater, MH= matric head, K= conductivity and FLW= chosen flow. In an iteration loop the correct flow is found. The integration goes at most over four intervals: [0,45], [45,170], [170,330] and [330,MH-rootzone] (last one on logarithmic scale).

Author
C. Rappoldt M. Wopereis
Date
January 1986, revised June 1990

@documentation Chapter 15 in documentation WOFOST Version 4.1 (1988)

281{
282 //calculation of pf from matrix head
283 double pf( std::log10( _cp));
284
285 //in case of small matric head (high water contents)
286 if ( cbm::flt_greater( 1.0, pf))
287 {
288 _flow = 0.0;
289 return LDNDC_ERR_OK;
290 }
291
292 double elog10( 2.302585);
293 double logst4( 2.518514);
294
295 double start[4] = {0.0, 45.0, 170.0, 330.0};
296
297 //number and width of integration intervals
298 int iint( 0);
299 double del[4] = { 0.0, 0.0, 0.0, 0.0};
300 for (int i1 = 0; i1 < 4; i1++)
301 {
302 if (i1 < 3)
303 {
304 del[i1] = std::min( start[i1+1], _cp) - start[i1];
305 }
306 else
307 {
308 del[i1] = pf - logst4;
309 }
310 if (del[i1] <= 0.0)
311 {
312 break;
313 }
314 iint += 1;
315 }
316
317 //preparation of three-point gaussian integration
318 double hulp[12] = { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
319 double conduc[12] = { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
320 for (int i1 = 0; i1 < iint; i1++)
321 {
322 for (int i2 = 0; i2 < 3; i2++)
323 {
324 double pgau[3] = {0.1127016654, 0.5, 0.8872983346};
325 double wgau[3] = {0.2777778, 0.4444444, 0.2777778};
326
327 int i3( 3 * i1 + i2);
328
329 //the three points in the full-width intervals are standard
330 double pfstan[9] = {0.705143, 1.352183, 1.601282, 1.771497, 2.031409, 2.192880, 2.274233, 2.397940, 2.494110};
331 double pfgau( pfstan[i3]);
332
333 //the three points in the last interval are calculated
334 if ( i1 == (iint-1))
335 {
336 if ( iint < 4)
337 {
338 pfgau = std::log10( start[iint-1] + pgau[i2] * del[iint-1]);
339 }
340 else
341 {
342 pfgau = logst4 + pgau[i2] * del[iint-1];
343 }
344 }
345
346 double const wcl_loc( ldndc::hydrology::water_content( std::exp(elog10 * pfgau), wc_.vga_sl[_sl], wc_.vgn_sl[_sl], wc_.vgm_sl[_sl], _saturation, 0.1));
347 conduc[i3] = ldndc::hydraulic_conductivity( wcl_loc, wc_.vgm_sl[_sl], wc_.mvg_tau_sl[_sl],
348 _saturation, wc_.wc_wp_sl[_sl], kst_sl[_sl] * cbm::CM_IN_M);
349 hulp[i3] = del[i1] * wgau[i2] * conduc[i3];
350 if (i3 > 8)
351 {
352 hulp[i3] = hulp[i3] * elog10 * std::exp( elog10 * pfgau);
353 }
354 }
355 }
356
357 //setting upper and lower limit
358 double const wcl_loc( ldndc::hydrology::water_content( _cp, wc_.vga_sl[_sl], wc_.vgn_sl[_sl], wc_.vgm_sl[_sl], _saturation, 0.1));
359 double const kms( ldndc::hydraulic_conductivity( wcl_loc, wc_.vgm_sl[_sl], wc_.mvg_tau_sl[_sl],
360 _saturation, wc_.wc_wp_sl[_sl], kst_sl[_sl] * cbm::CM_IN_M));
361 double fu( (_cp <= _delta_z) ? 0.0 : 1.27);
362 double fl( (_cp >= _delta_z) ? 0.0 : -kms);
363
364 if ( !cbm::flt_equal( _cp, _delta_z))
365 {
366 //iteration loop
367 int imax( 3 * iint);
368 for ( int i1 = 0; i1 < 14; i1++)
369 {
370 double flw( (fu + fl) / 2.0);
371 double df( (fu - fl) / 2.0);
372 if ( (df < 0.01) &&
373 ((df / std::abs(flw)) < 0.1))
374 {
375 break;
376 }
377
378 double z( 0.0);
379 for ( int i2 = 0; i2 < imax; i2++)
380 {
381 z = z + hulp[i2] / (conduc[i2] + flw);
382 }
383 if (z >= _delta_z)
384 {
385 fl = flw;
386 }
387 if (z <= _delta_z)
388 {
389 fu = flw;
390 }
391 }
392 }
393
394 //flow output in mm
395 _flow = 10.0 * (fu + fl) / 2.0;
396
397 return LDNDC_ERR_OK;
398}
ldndc::EcHy::kst_sl
lvector< double > kst_sl
Definition echy.h:113

References kst_sl.

Referenced by EcHyEvapotranspiration().

Here is the caller graph for this function:

◆ EcHyFlood()

lerr_t ldndc::EcHy::EcHyFlood ( )
private

sets hydrologic conditions during flooding events, e.g.,

Flooding.

  • surface water table
  • bund height
18{
19 lerr_t rc = m_eventflood.solve();
20 if ( rc != LDNDC_ERR_OK){ return rc; }
21
22 bund_height = m_eventflood.get_bund_height();
23
24 /* todo */
25 //irrigation_amount = m_eventflood.get_irrigation_amount();
26
27 /* todo */
28 //double const irrigation_height = m_eventflood.get_irrigation_height();
29 double const max_percolation = m_eventflood.get_maximum_percolation();
30 double const water_table_flooding = m_eventflood.get_water_table();
31
32 irrigation_switch = NONE;
33 minimum_watertable_height = 0.0;
34 if ( cbm::is_valid( water_table_flooding))
35 {
36 if ( cbm::flt_greater_zero( water_table_flooding))
37 {
38 irrigation_switch = CONSTANT_POSITIVE_WATER_TABLE;
39 minimum_watertable_height = water_table_flooding;
40 bund_height = water_table_flooding;
41 }
42 else
43 {
44 irrigation_switch = CONSTANT_NEGATIVE_WATER_TABLE;
45 minimum_watertable_height = water_table_flooding;
46 bund_height = 0.0;
47 }
48 }
49
50 bool cracked( false);
51 for (size_t sl = 0; sl < soillayers_in->soil_layer_cnt()-1; sl++)
52 {
53 if ( cbm::flt_greater_equal( sc_.depth_sl[sl], 0.3))
54 {
55 if ( cbm::is_valid( max_percolation))
56 {
57 if ((kst_eq_sl[sl] > max_percolation / lclock()->time_resolution()))
58 {
59
60 kst_sl[sl] = cbm::bound_min( 0.0,
61 kst_eq_sl[sl] -
62 (kst_eq_sl[sl] - max_percolation / lclock()->time_resolution()) *
63 (sc_.depth_sl[sl] - 0.2) / (sc_.depth_sl[soillayers_in->soil_layer_cnt()-1] - 0.2));
64 }
65 kst_sl[sl] = max_percolation / lclock()->time_resolution();
66 }
67 else if ( wl_sl[sl] < -0.8 * pore_space( sl))
68 {
69 cracked = true;
70 }
71
72 if ( cracked)
73 {
74 double const time_rate( 0.1 / lclock()->time_resolution());
75 kst_sl[sl] = cbm::bound_min( 0.0,
76 kst_sl[sl] - (kst_sl[sl] - kst_eq_sl[sl]) * time_rate);
77 }
78 break;
79 }
80 }
81
82 return LDNDC_ERR_OK;
83}

References kst_sl, and wl_sl.

Referenced by solve().

Here is the caller graph for this function:

◆ EcHyGroundwater()

lerr_t ldndc::EcHy::EcHyGroundwater ( )
private

...

Negative groundwater table represents water on the soil surface

20{
21 cbm::invalidate( kst_bottom_gw_vertical);
22 if ( cbm::flt_greater( sc_.depth_sl[soillayers_in->soil_layer_cnt()-1], gw_depth_static))
23 {
24 // negative groundwater table represents water above soil surface
25 double const surface_water_new( surface_water + throughfall_water + surface_snowmelt + irrigation
26 - ev_surfacewater - wlfl_sl[0] - runoff - bypass_fill_sl.sum());
27 if ( cbm::flt_greater_zero( - gw_depth_static))
28 {
29 gw_fill_surface = - gw_depth_static - surface_water_new;
30 }
31 // remove surface water if static groundwater existed in previous time step
32 else if ( cbm::flt_greater_zero( - gw_depth_static_old))
33 {
34 gw_fill_surface = - cbm::bound_max( surface_water, surface_water_new);
35 }
36
37 //groundwater access
38 {
39 size_t sl( 0);
40 double layer_midpoint_depth( 0.5 * sc_.h_sl[sl]);
41 if ( cbm::flt_greater_equal( layer_midpoint_depth, gw_depth_static))
42 {
43 gw_fill_sl[sl] = cbm::bound_min( 0.0,
44 0.999 * pore_space( sl) - wc_.ice_sl[sl] * sc_.h_sl[sl]
45 - (wl_sl[sl] + EcHySoilWaterChange( sl)));
46 wc_.accumulated_groundwater_access += gw_fill_sl[sl];
47 }
48 else
49 {
50 double const wl_new( wl_sl[sl] + EcHySoilWaterChange( sl));
51 double const wc_new( cbm::bound( wc_.wc_res_sl[sl], wl_new / sc_.h_sl[sl], wc_.wc_sat_sl[sl]));
52 kst_bottom_gw_vertical = cbm::bound_min( 0.0,
53 ldndc::hydraulic_conductivity(
54 wc_new,
55 wc_.vgm_sl[sl],
56 wc_.mvg_tau_sl[sl],
57 wc_.wc_sat_sl[sl],
58 wc_.wc_res_sl[sl],
59 kst_sl[sl]));
60 gw_fill_sl[sl] = 0.0;
61 }
62
63 for ( sl = 1; sl < soillayers_in->soil_layer_cnt()-1; sl++)
64 {
65 layer_midpoint_depth = sc_.depth_sl[sl] - 0.5 * sc_.h_sl[sl];
66 if ( cbm::flt_greater_equal( layer_midpoint_depth, gw_depth_static))
67 {
68 gw_fill_sl[sl] = cbm::bound_min( 0.0,
69 0.999 * pore_space( sl) - wc_.ice_sl[sl] * sc_.h_sl[sl]
70 - (wl_sl[sl] + EcHySoilWaterChange( sl)));
71 wc_.accumulated_groundwater_access += gw_fill_sl[sl];
72 }
73 else
74 {
75 double const wl_new( wl_sl[sl] + EcHySoilWaterChange( sl));
76 double const wc_new( cbm::bound( wc_.wc_res_sl[sl], wl_new / sc_.h_sl[sl], wc_.wc_sat_sl[sl]));
77 kst_bottom_gw_vertical = cbm::bound_min( 0.0,
78 ldndc::hydraulic_conductivity(
79 wc_new,
80 wc_.vgm_sl[sl],
81 wc_.mvg_tau_sl[sl],
82 wc_.wc_sat_sl[sl],
83 wc_.wc_res_sl[sl],
84 kst_sl[sl]));
85 gw_fill_sl[sl] = 0.0;
86 }
87 }
88
89 sl = soillayers_in->soil_layer_cnt()-1;
90 layer_midpoint_depth = sc_.depth_sl[sl] - 0.5 * sc_.h_sl[sl];
91 if ( cbm::flt_greater_equal( layer_midpoint_depth, gw_depth_static))
92 {
93 gw_fill_sl[sl] = cbm::bound_min( 0.0,
94 0.999 * pore_space( sl) - wc_.ice_sl[sl] * sc_.h_sl[sl]
95 - (wl_sl[sl] + EcHySoilWaterChange( sl)));
96 wc_.accumulated_groundwater_access += gw_fill_sl[sl];
97 }
98 else
99 {
100 double const wl_new( wl_sl[sl] + EcHySoilWaterChange( sl));
101 double const wc_new( cbm::bound( wc_.wc_res_sl[sl], wl_new / sc_.h_sl[sl], wc_.wc_sat_sl[sl]));
102 kst_bottom_gw_vertical = cbm::bound_min( 0.0,
103 ldndc::hydraulic_conductivity(
104 wc_new,
105 wc_.vgm_sl[sl],
106 wc_.mvg_tau_sl[sl],
107 wc_.wc_sat_sl[sl],
108 wc_.wc_res_sl[sl],
109 kst_sl[sl]));
110 gw_fill_sl[sl] = 0.0;
111 }
112 }
113 }
114
115 //groundwater removal
116 for (size_t sl = 1; sl < soillayers_in->soil_layer_cnt()-1; sl++)
117 {
118 double const layer_midpoint_depth( sc_.depth_sl[sl]);
119 if ( cbm::flt_greater_equal( layer_midpoint_depth, gw_depth_static_old) &&
120 cbm::flt_less( layer_midpoint_depth, gw_depth_static) &&
121 cbm::flt_greater( wl_sl[sl] + EcHySoilWaterChange( sl), wc_.wc_fc_sl[sl] * sc_.h_sl[sl]))
122 {
123 double wl_reduce( wl_sl[sl] + EcHySoilWaterChange( sl) - wc_.wc_fc_sl[sl] * sc_.h_sl[sl]);
124 gw_fill_sl[sl] -= wl_reduce;
125 }
126 }
127
128 //capillary rise
129 if ( have_capillary_flow)
130 {
131 size_t const STEPS( 10);
132 for (size_t steps = 0; steps < STEPS; steps++)
133 {
134 for (size_t sl = 0; sl < soillayers_in->soil_layer_cnt()-1; sl++)
135 {
136 double const layer_midpoint_depth( sc_.depth_sl[sl] - 0.5 * sc_.h_sl[sl]);
137 if ( cbm::flt_greater_equal( layer_midpoint_depth, gw_depth_static))
138 {
139 break;
140 }
141
142 /* Water contents with until now calculated water balance */
143 double const wl_new_0( wl_sl[sl] + EcHySoilWaterChange( sl));
144 double const wl_new_1( wl_sl[sl+1] + EcHySoilWaterChange( sl+1));
145
146 double const wc_vg_0( cbm::bound( wc_.wc_res_sl[sl], wl_new_0 / sc_.h_sl[sl], wc_.wc_sat_sl[sl]));
147 double const wc_vg_1( cbm::bound( wc_.wc_res_sl[sl+1], wl_new_1 / sc_.h_sl[sl+1], wc_.wc_sat_sl[sl+1]));
148
149 /* Unsaturated hydraulic conductivity with until now calculated water balance */
150 double const kust_0( cbm::bound_min( 0.0,
151 ldndc::hydraulic_conductivity(
152 wc_vg_0, wc_.vgm_sl[sl], wc_.mvg_tau_sl[sl],
153 wc_.wc_sat_sl[sl], wc_.wc_res_sl[sl],
154 1.0/double(STEPS)*kst_sl[sl])));
155 double const kust_1( cbm::bound_min( 0.0,
156 ldndc::hydraulic_conductivity(
157 wc_vg_1, wc_.vgm_sl[sl+1], wc_.mvg_tau_sl[sl],
158 wc_.wc_sat_sl[sl+1], wc_.wc_res_sl[sl+1],
159 1.0/double(STEPS)*kst_sl[sl+1])));
160
161 double const capillary_pressure_0( ldndc::capillary_pressure(
162 wl_new_0/sc_.h_sl[sl], wc_.vga_sl[sl],
163 wc_.vgn_sl[sl], wc_.vgm_sl[sl],
164 wc_.wc_sat_sl[sl], wc_.wc_res_sl[sl]));
165 double const capillary_pressure_1( ldndc::capillary_pressure(
166 wl_new_1/sc_.h_sl[sl+1], wc_.vga_sl[sl+1],
167 wc_.vgn_sl[sl+1], wc_.vgm_sl[sl+1],
168 wc_.wc_sat_sl[sl], wc_.wc_res_sl[sl+1]));
169
170 double const kust( cbm::flt_greater( capillary_pressure_0, capillary_pressure_1) ? kust_1 : kust_0 );
171 double const flow( (capillary_pressure_1 - capillary_pressure_0) / (0.5 * (sc_.h_sl[sl] + sc_.h_sl[sl+1])) * kust);
172
173 /* downward flow*/
174 if ( cbm::flt_greater_zero( flow))
175 {
176 double const left( wl_sl[sl] + EcHySoilWaterChange( sl) - (wc_.wc_res_sl[sl] * sc_.h_sl[sl]));
177 double const space( 0.999 * pore_space( sl+1) - wc_.ice_sl[sl+1] * sc_.h_sl[sl+1]
178 - (wl_sl[sl+1] + EcHySoilWaterChange( sl+1)));
179
180 if ( cbm::flt_greater_zero( left) &&
181 cbm::flt_greater_zero( space))
182 {
183 double const bound_flow( cbm::bound_max( flow, std::min( left, space)) / static_cast<double>(STEPS));
184 cr_fill_sl[sl] -= bound_flow;
185 cr_fill_sl[sl+1] += bound_flow;
186 cr_down_sl[sl] += bound_flow;
187 }
188 }
189 /* upward flow */
190 else if ( cbm::flt_less( flow, 0.0))
191 {
192 double const left( wl_sl[sl+1] + EcHySoilWaterChange( sl+1) - (wc_.wc_res_sl[sl+1] * sc_.h_sl[sl+1]));
193 double const space( 0.999 * pore_space( sl) - wc_.ice_sl[sl] * sc_.h_sl[sl]
194 - (wl_sl[sl] + EcHySoilWaterChange( sl)));
195
196 if ( cbm::flt_greater_zero( left) &&
197 cbm::flt_greater_zero( space))
198 {
199 double const bound_flow( cbm::bound_max( -flow, std::min( left, space)) / static_cast<double>(STEPS));
200 cr_fill_sl[sl] += bound_flow;
201 cr_fill_sl[sl+1] -= bound_flow;
202 cr_up_sl[sl] += bound_flow;
203 }
204 }
205 }
206
207 size_t sl( soillayers_in->soil_layer_cnt()-1);
208 double const layer_midpoint_depth( sc_.depth_sl[sl] - 0.5 * sc_.h_sl[sl]);
209 if ( have_capillary_rise &&
210 //groundwater table below discretized soil
211 cbm::flt_greater( gw_depth_static, layer_midpoint_depth) &&
212 //groundwater not more than 10m below discretized soil
213 cbm::flt_less( gw_depth_static, layer_midpoint_depth + 10))
214 {
215 /* Water contents with until now calculated water balance */
216 double const wl_new( wl_sl[sl] + EcHySoilWaterChange( sl));
217
218 double const wc_vg( cbm::bound( wc_.wc_res_sl[sl], wl_new / sc_.h_sl[sl], wc_.wc_sat_sl[sl]));
219
220 /* Unsaturated hydraulic conductivity with until now calculated water balance */
221 double const kust( cbm::bound_min( 0.0,
222 ldndc::hydraulic_conductivity(
223 wc_vg, wc_.vgm_sl[sl], wc_.mvg_tau_sl[sl],
224 wc_.wc_sat_sl[sl], wc_.wc_res_sl[sl],
225 1.0/double(STEPS)*kst_sl[sl])));
226
227 double const capillary_pressure( ldndc::capillary_pressure(
228 wl_new/sc_.h_sl[sl], wc_.vga_sl[sl],
229 wc_.vgn_sl[sl], wc_.vgm_sl[sl],
230 wc_.wc_sat_sl[sl], wc_.wc_res_sl[sl]));
231
232 /* upward flow */
233 double const flow( -capillary_pressure / (gw_depth_static - layer_midpoint_depth) * kust);
234 if ( cbm::flt_less( flow, 0.0))
235 {
236 double const space( 0.999 * pore_space( sl) - wc_.ice_sl[sl] * sc_.h_sl[sl]
237 - (wl_sl[sl] + EcHySoilWaterChange( sl)));
238 if ( cbm::flt_greater_zero( space))
239 {
240 double const bound_flow( cbm::bound_max( -flow, 1.0 / static_cast<double>(STEPS) * space));
241 cr_fill_sl[sl] += bound_flow;
242 cr_fill_groundwater += bound_flow;
243 cr_up_sl[sl] += bound_flow;
244 }
245 }
246 }
247 }
248 }
249
250 return LDNDC_ERR_OK;
251}
ldndc::EcHy::bypass_fill_sl
lvector< double > bypass_fill_sl
Definition echy.h:104
ldndc::EcHy::EcHySoilWaterChange
double EcHySoilWaterChange(size_t)
Definition echy.cpp:511
ldndc::EcHy::cr_fill_sl
lvector< double > cr_fill_sl
Definition echy.h:96
ldndc::EcHy::gw_depth_static
double gw_depth_static
Definition echy.h:74

References bypass_fill_sl, cr_fill_groundwater, cr_fill_sl, EcHySoilWaterChange(), ev_surfacewater, gw_depth_static, gw_fill_sl, gw_fill_surface, kst_sl, runoff, throughfall_water, wl_sl, and wlfl_sl.

Referenced by EcHyEvapotranspiration(), and solve().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ EcHyIntegration()

lerr_t ldndc::EcHy::EcHyIntegration ( )
private

...

Integrates all water state variables.

524{
525 canopy_snow = cbm::bound_min( 0.0,
526 canopy_snow - canopy_snowmelt);
527 surface_snow = cbm::bound_min( 0.0,
528 surface_snow
529 + throughfall_snow
530 - surface_snowmelt);
531 if ( cbm::flt_less( surface_snow, 1.0e-9))
532 {
533 surface_snow = 0.0;
534 }
535 canopy_water = cbm::bound_min( 0.0,
536 canopy_water
537 + canopy_snowmelt);
538 if ( cbm::flt_less( canopy_water, 1.0e-9))
539 {
540 canopy_water = 0.0;
541 }
542 surface_water = cbm::bound_min( 0.0,
543 surface_water
544 + throughfall_water + surface_snowmelt + gw_fill_surface + irrigation
545 - ev_surfacewater - wlfl_sl[0] - runoff - bypass_fill_sl.sum());
546 if ( cbm::flt_less( surface_water, 1.0e-9))
547 {
548 surface_water = 0.0;
549 }
550
551 for (size_t sl = 0; sl < soillayers_in->soil_layer_cnt(); sl++)
552 {
553 wl_sl[sl] = cbm::bound_min( 0.0, wl_sl[sl] + EcHySoilWaterChange( sl));
554 wc_.mskpa_sl[sl] = ldndc::hydrology::capillary_pressure( wl_sl[sl]/sc_.h_sl[sl], wc_.vga_sl[sl],
555 wc_.vgn_sl[sl], wc_.vgm_sl[sl],
556 wc_.wc_sat_sl[sl], wc_.wc_res_sl[sl]);
557 }
558
559 return LDNDC_ERR_OK;
560}
ldndc::EcHy::canopy_snowmelt
double canopy_snowmelt
Definition echy.h:144

References bypass_fill_sl, canopy_snowmelt, canopy_water, EcHySoilWaterChange(), ev_surfacewater, gw_fill_surface, runoff, throughfall_snow, throughfall_water, wl_sl, and wlfl_sl.

Referenced by EcHyEvapotranspiration(), and solve().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ EcHySoilWaterChange()

double ldndc::EcHy::EcHySoilWaterChange ( size_t _sl)
private

Integrates all water state variables.

513{
514 return wlfl_sl[_sl] - wlfl_sl[_sl+1] + cr_fill_sl[_sl] - trwl_sl[_sl] - evsws_sl[_sl] + gw_fill_sl[_sl];
515}

References cr_fill_sl, evsws_sl, gw_fill_sl, trwl_sl, and wlfl_sl.

Referenced by EcHyEvapotranspiration(), EcHyGetAvailableWaterEvaporation(), EcHyGetAvailableWaterTranspiration(), EcHyGroundwater(), EcHyIntegration(), and EcHyIrrigation().

Here is the caller graph for this function:

◆ solve()

lerr_t ldndc::EcHy::solve ( )

Kicks off computation for one time step.

40{
41 /* reset of all fluxes */
42 EcHyreset();
43
44 /* flooding events */
45 lerr_t rc = EcHyFlood();
46 if ( rc != LDNDC_ERR_OK)
47 {
48 KLOGERROR("Handling flooding event failed in: ", name());
49 return rc;
50 }
51
52 /* tasks:
53 * - update internal state (e.g., throughfall_water, ...)
54 * - calculate evapotranspiration on demand
55 */
56 rc = EcHyStepInit();
57 if ( rc != LDNDC_ERR_OK)
58 { return rc; }
59
60 /* begin water balance */
61 double balance( -1.0);
62 EcHyBalanceCheck( balance);
63
64 /* snow and ice calculated */
65 EcHySnowIce();
66
67 /*
68 * apply irrigation triggered by
69 * - irrigation events
70 * - flooding events
71 */
72 EcHyIrrigation();
73
74 double const potential_infiltration( cbm::flt_greater( sc_.depth_sl[0], gw_depth_static) ?
75 irrigation + throughfall_water + surface_snowmelt:
76 irrigation + throughfall_water + surface_snowmelt + surface_water);
77
78 if ( cbm::flt_greater_zero( potential_infiltration))
79 {
80 /* actual/potential runoff */
81 double potential_runoff( 0.0);
82 if ( runoff_method == "curvenumber")
83 {
84 if ( cbm::flt_greater_zero( bund_height) &&
85 cbm::flt_greater( potential_infiltration, bund_height))
86 {
87 runoff = EcHyRunoff( potential_infiltration - bund_height);
88 }
89 else
90 {
91 runoff = EcHyRunoff( potential_infiltration);
92 }
93 }
94 else
95 {
96 if ( cbm::flt_greater_zero( bund_height) &&
97 cbm::flt_greater( potential_infiltration, bund_height))
98 {
99 double const f_runoff( cbm::bound_max( m_param->FRUNOFF() / lclock()->time_resolution(), 1.0));
100 potential_runoff = cbm::bound_min( 0.0, f_runoff * (potential_infiltration - bund_height));
101 }
102 else
103 {
104 double const f_runoff( cbm::bound_max( m_param->FRUNOFF() / lclock()->time_resolution(), 1.0));
105 potential_runoff = cbm::bound_min( 0.0, f_runoff * potential_infiltration);
106 }
107 }
108
109 /* potential bypass */
110 double const potential_bypass_flow( m_param->CRACK_FRACTION() * potential_infiltration);
111
112 double const potential_flows_sum( potential_runoff + potential_bypass_flow + kst_sl[0]);
113
114 if ( cbm::flt_greater_zero( potential_flows_sum))
115 {
116 /* runoff */
117 double max_flow( cbm::bound_max( potential_runoff, potential_infiltration));
118 runoff += cbm::bound_max( potential_runoff / potential_flows_sum * potential_infiltration,
119 max_flow);
120
121 /* bypass */
122 max_flow = cbm::bound_max( potential_bypass_flow, potential_infiltration - runoff);
123 EcHyBypassFlow( m_param->CRACK_DEPTH(), cbm::bound_max( potential_bypass_flow / potential_flows_sum * potential_infiltration,
124 max_flow));
125
126 /* infiltration */
127 wlfl_sl[0] = cbm::bound_max( potential_infiltration - runoff - bypass_fill_sl.sum(),
128 kst_sl[0]);
129 }
130 }
131
132 for ( size_t sl = 0; sl < soillayers_in->soil_layer_cnt(); ++sl)
133 {
134 EcHyPercolation( 1, sl);
135 }
136
137 if ( cbm::is_valid( kst_bottom))
138 {
139 wlfl_sl[soillayers_in->soil_layer_cnt()] = cbm::bound_max( wlfl_sl[soillayers_in->soil_layer_cnt()],
140 kst_bottom);
141 }
142 if ( cbm::is_valid( kst_bottom_gw_lateral))
143 {
144 wlfl_sl[soillayers_in->soil_layer_cnt()] = cbm::bound_max( wlfl_sl[soillayers_in->soil_layer_cnt()],
145 kst_bottom_gw_lateral);
146 }
147 if ( cbm::is_valid( kst_bottom_gw_vertical))
148 {
149 wlfl_sl[soillayers_in->soil_layer_cnt()] = cbm::bound_max( wlfl_sl[soillayers_in->soil_layer_cnt()],
150 kst_bottom_gw_vertical);
151 }
152
153 if ( cbm::flt_greater( sc_.depth_sl[soillayers_in->soil_layer_cnt()-1], gw_depth_static-1.0))
154 {
155 if ( groundwater_percolation_method == "relative")
156 {
157 wlfl_sl[soillayers_in->soil_layer_cnt()] = cbm::bound_max( wlfl_sl[soillayers_in->soil_layer_cnt()],
158 m_param->GROUNDWATER_PERCOLATION_RELATIVE() *
159 kst_sl[soillayers_in->soil_layer_cnt()-1]);
160 }
161 else if ( groundwater_percolation_method == "absolute")
162 {
163 wlfl_sl[soillayers_in->soil_layer_cnt()] = cbm::bound_max( wlfl_sl[soillayers_in->soil_layer_cnt()],
164 m_param->GROUNDWATER_PERCOLATION_ABSOLUTE() /
165 lclock()->time_resolution() * cbm::M_IN_MM);
166 }
167 else
168 {
169 wlfl_sl[soillayers_in->soil_layer_cnt()] = cbm::bound_max( wlfl_sl[soillayers_in->soil_layer_cnt()],
170 gw_recharge_mean);
171 }
172 }
173
174 EcHyEvapotranspiration();
175
176 /* upward water flow */
177 for ( size_t sl = soillayers_in->soil_layer_cnt(); sl-- > 0; )
178 {
179 EcHyPercolation( 2, sl);
180 }
181
182 /* groundwater interaction:
183 * - set soil layers water-saturated
184 * - capillary rise
185 */
186 EcHyGroundwater();
187
188 /* integration of state variables */
189 EcHyIntegration();
190
191 /* update external state */
192 EcHyStepExit();
193
194 /* perform water balance */
195 rc = EcHyBalanceCheck( balance);
196 if ( rc){ return rc; }
197
198 rc = EcHy_check_for_negative_value("exit");
199 if ( rc){ return rc; }
200
201 return LDNDC_ERR_OK;
202}
ldndc::EcHy::EcHyPercolation
lerr_t EcHyPercolation(size_t, size_t)
Calculates water percolation within the soil profile.
ldndc::EcHy::EcHySnowIce
lerr_t EcHySnowIce()
Calls SnowDNDC for the calculation of snowpack and soil ice formation.
Definition echy-snow.cpp:21
ldndc::EcHy::EcHyFlood
lerr_t EcHyFlood()
sets hydrologic conditions during flooding events, e.g.,
Definition echy-management.cpp:17
ldndc::EcHy::EcHyStepInit
lerr_t EcHyStepInit()
...
Definition echy.cpp:244
ldndc::EcHy::EcHyEvapotranspiration
lerr_t EcHyEvapotranspiration()
Calculates evapotranspiration within the soil profile.
ldndc::EcHy::EcHyGroundwater
lerr_t EcHyGroundwater()
...
Definition echy-groundwater.cpp:19
ldndc::EcHy::EcHyBalanceCheck
lerr_t EcHyBalanceCheck(double &)
...
Definition echy.cpp:568
ldndc::EcHy::EcHyIrrigation
lerr_t EcHyIrrigation()
Irrigation.
Definition echy-management.cpp:91
ldndc::EcHy::EcHyreset
void EcHyreset()
...
Definition echy.cpp:210
ldndc::EcHy::EcHyBypassFlow
lerr_t EcHyBypassFlow(double, double)
Calculates water percolation within the soil profile.
ldndc::EcHy::EcHyIntegration
lerr_t EcHyIntegration()
...
Definition echy.cpp:523
ldndc::EcHy::EcHyStepExit
lerr_t EcHyStepExit()
...
Definition echy.cpp:383
ldndc::EcHy::kst_bottom
double kst_bottom
Definition echy.h:69
ldndc::EcHy::runoff_method
cbm::string_t runoff_method
Definition echy.h:200

References bypass_fill_sl, EcHyBalanceCheck(), EcHyBypassFlow(), EcHyEvapotranspiration(), EcHyFlood(), EcHyGroundwater(), EcHyIntegration(), EcHyIrrigation(), EcHyPercolation(), EcHyreset(), EcHySnowIce(), EcHyStepExit(), EcHyStepInit(), gw_depth_static, kst_bottom, kst_sl, runoff, runoff_method, throughfall_water, and wlfl_sl.

Here is the call graph for this function:

Member Data Documentation

◆ accumulated_irrigation_old

double ldndc::EcHy::accumulated_irrigation_old
private

Accumulated irrigation of last time step [m]

Referenced by EcHyIrrigation().

◆ accumulated_potentialtranspiration_old

double ldndc::EcHy::accumulated_potentialtranspiration_old
private

Accumulated potential transpiration of last time step [m]

Referenced by EcHyPotentialEvapotranspiration().

◆ bypass_fill_sl

lvector< double > ldndc::EcHy::bypass_fill_sl
private

Soillayer water input due to bypass flow

Referenced by EcHyGroundwater(), EcHyIntegration(), EcHyreset(), and solve().

◆ canopy_snowmelt

double ldndc::EcHy::canopy_snowmelt
private

Snow melting

Referenced by EcHyIntegration(), and EcHySnowIce().

◆ canopy_water

double ldndc::EcHy::canopy_water
private

Amount of canopy intercepted water

Referenced by EcHyBalanceCheck(), EcHyCalculateLeafWaterDistribution(), EcHyIntegration(), EcHyPotentialEvapotranspiration(), EcHySnowIce(), and EcHyStepInit().

◆ cr_fill_groundwater

double ldndc::EcHy::cr_fill_groundwater
private

Water addition to last soil layer from groundwater water by cappillary rise (only used for water balance calculation)

Referenced by EcHyBalanceCheck(), EcHyGroundwater(), and EcHyreset().

◆ cr_fill_sl

lvector< double > ldndc::EcHy::cr_fill_sl
private

Soillayer water input due to capillary rise

Referenced by EcHyGroundwater(), EcHyreset(), and EcHySoilWaterChange().

◆ daily_potential_evapotranspiration

double ldndc::EcHy::daily_potential_evapotranspiration
private

...

Referenced by EcHyPotentialEvapotranspiration(), and EcHyStepExit().

◆ daily_potential_leaf_evaporation

double ldndc::EcHy::daily_potential_leaf_evaporation
private

...

Referenced by EcHyPotentialEvapotranspiration().

◆ daily_potential_soil_evaporation

double ldndc::EcHy::daily_potential_soil_evaporation
private

...

Referenced by EcHyPotentialEvapotranspiration().

◆ daily_potential_transpiration

double ldndc::EcHy::daily_potential_transpiration
private

...

Referenced by EcHyPotentialEvapotranspiration().

◆ ev_leaf

double ldndc::EcHy::ev_leaf
private

Evaporation from plant surface

Referenced by EcHyBalanceCheck(), and EcHyStepExit().

◆ ev_surfacewater

double ldndc::EcHy::ev_surfacewater
private

Evaporation from surface water

Referenced by EcHyBalanceCheck(), EcHyGroundwater(), EcHyIntegration(), EcHyreset(), and EcHyStepExit().

◆ evapotranspiration_method

cbm::string_t ldndc::EcHy::evapotranspiration_method
private

...

Referenced by EcHyPotentialEvapotranspiration().

◆ evsws_sl

lvector< double > ldndc::EcHy::evsws_sl
private

Evaporation from soil layer

Referenced by EcHyBalanceCheck(), EcHyreset(), EcHySoilWaterChange(), and EcHyStepExit().

◆ gw_depth_static

double ldndc::EcHy::gw_depth_static
private

Depth of groundwater table [m]

Referenced by EcHyGroundwater(), EcHyStepExit(), EcHyStepInit(), and solve().

◆ gw_fill_sl

lvector< double > ldndc::EcHy::gw_fill_sl
private

Soillayer water input due to groundwater flow

Referenced by EcHyBalanceCheck(), EcHyGroundwater(), EcHyreset(), and EcHySoilWaterChange().

◆ gw_fill_surface

double ldndc::EcHy::gw_fill_surface
private

Water addition to surface water from groundwater water boundary condition

Referenced by EcHyBalanceCheck(), EcHyGroundwater(), EcHyIntegration(), and EcHyreset().

◆ kst_bottom

double ldndc::EcHy::kst_bottom
private

Saturated hydraulic conductivity below last soil layer [cm:min-1]

Referenced by EcHyStepInit(), and solve().

◆ kst_sl

lvector< double > ldndc::EcHy::kst_sl
private

Saturated hydraulic conductivity

Referenced by EcHyCapillaryRise(), EcHyFlood(), EcHyGroundwater(), and solve().

◆ m_icecontent_in

WaterCycleSnowDNDC::IceContentStateIn ldndc::EcHy::m_icecontent_in
private

...

Referenced by EcHySnowIce().

◆ runoff

double ldndc::EcHy::runoff
private

Lateral surface runoff

Referenced by EcHyBalanceCheck(), EcHyGroundwater(), EcHyIntegration(), EcHyreset(), EcHyStepExit(), and solve().

◆ runoff_method

cbm::string_t ldndc::EcHy::runoff_method
private

...

Referenced by solve().

◆ thornthwaite_heat_index

double ldndc::EcHy::thornthwaite_heat_index
private

...

Referenced by EcHyPotentialEvapotranspiration().

◆ throughfall_snow

double ldndc::EcHy::throughfall_snow
private

Snowfall

Referenced by EcHyBalanceCheck(), EcHyIntegration(), EcHyreset(), and EcHySnowIce().

◆ throughfall_water

double ldndc::EcHy::throughfall_water
private

Throughfall to soil

Referenced by EcHyBalanceCheck(), EcHyGroundwater(), EcHyIntegration(), EcHySnowIce(), EcHyStepExit(), EcHyStepInit(), and solve().

◆ trwl_sl

lvector< double > ldndc::EcHy::trwl_sl
private

Water withdrawl by transpiration [m d−1]

Referenced by EcHyBalanceCheck(), EcHyreset(), EcHySoilWaterChange(), and EcHyStepExit().

◆ wl_sl

lvector< double > ldndc::EcHy::wl_sl
private

Amount of water in soil layer

Referenced by EcHyBalanceCheck(), EcHyFlood(), EcHyGetAvailableWaterEvaporation(), EcHyGetAvailableWaterTranspiration(), EcHyGroundwater(), EcHyIntegration(), EcHyIrrigation(), EcHySnowIce(), EcHyStepExit(), and EcHyStepInit().

◆ wlfc_sl

lvector< double > ldndc::EcHy::wlfc_sl
private

Amount of water in soil layer at field capacity

Referenced by EcHyIrrigation().

◆ wlfl_sl

lvector< double > ldndc::EcHy::wlfl_sl
private

Water flux at boundaries of soil layer [m d-1] the flux of water between soil layer is tracked by wlfl (mm d−1). in total, there are nl+1 flow rates, where nl is the number of soil layers, wlfl[1] is the flow rate between the ponded water layer and the soil surface, wlfl[2] is the flow rate between soil layer 1 and soil layer 2, etc.

Referenced by EcHyBalanceCheck(), EcHyGroundwater(), EcHyIntegration(), EcHyreset(), EcHySoilWaterChange(), EcHyStepExit(), and solve().

◆ wlwp_sl

lvector< double > ldndc::EcHy::wlwp_sl
private

Amount of water in soil layer at wilting point

Referenced by EcHyGetWiltingPoint().