MARRMOT_44
Implementation
MARRMoT v2.1.1 rev eeb7e15 m_44_echo_16p_6s; standardized continuous structure, not original model
Time step: daily
Backend: octave-cli + Octave optim
Calibrated parameters: 16
Temperature required: yes
Parameters and initial configuration
Parameter |
Supported calibration range |
Default |
|---|---|---|
|
0 to 5 |
2.5 |
|
-3 to 5 |
1.0 |
|
-3 to 3 |
0.0 |
|
0 to 20 |
10.0 |
|
0 to 1 |
0.5 |
|
0 to 2 |
1.0 |
|
0 to 1 |
0.5 |
|
0 to 200 |
100.0 |
|
1 to 2000 |
1000.5 |
|
0.05 to 0.95 |
0.5 |
|
0.05 to 0.95 |
0.5 |
|
0 to 1 |
0.5 |
|
0 to 5 |
2.5 |
|
0 to 20 |
10.0 |
|
0 to 1 |
0.5 |
|
0 to 1 |
0.5 |
|
Fixed initial/configuration value |
0.0 |
|
Fixed initial/configuration value |
0.0 |
|
Fixed initial/configuration value |
0.0 |
|
Fixed initial/configuration value |
0.0 |
|
Fixed initial/configuration value |
0.0 |
|
Fixed initial/configuration value |
0.0 |
Ranges/defaults are implementation contracts, not universal priors or a recommended basin calibration. Consult source comments for parameter units and coupling.
Governing equations
The following source is the exact model kernel used by this adapter. For MARRMoT it includes the state derivative and each referenced flux function; the solver and routing are described above. Original notices and source citations are retained in the files.
function [dS, fluxes] = model_fun(obj, S)
% parameters
theta = obj.theta;
rho = theta(1); % Maximum interception storage [mm]
ts = theta(2); % Threshold temperature for snowfall [oC]
tm = theta(3); % Threshold temperature for snowmelt [oC]
as = theta(4); % Degree-day factor [mm/oC/d]
af = theta(5); % Refreezing reduction factor [-]
gmax = theta(6); % Maximum melt due to ground-heat flux [mm/d]
the = theta(7); % Water-holding capacity of snow [-]
phi = theta(8); % Maximum infiltration rate [mm/d]
smax = theta(9); % Maximum soil moisture storage [mm]
fsm = theta(10); % Plant stress point as a fraction of Smax [-]
fsw = theta(11); % Wilting point as fraction of sm [-]
ksat = theta(12); % Runoff rate from soil moisture [d-1]
c = theta(13); % Runoff non-linearity from soil moisture [-]
lmax = theta(14); % Groundwater flux [mm/d]
kf = theta(15); % Runoff coefficient [d-1]
ks = theta(16); % Runoff coefficient [d-1]
% auxiliary parameters
aux_theta = obj.aux_theta;
sm = aux_theta(1); % Plant stress point [mm]
sw = aux_theta(2); % Wilting point [mm]
% delta_t
delta_t = obj.delta_t;
% stores
S1 = S(1);
S2 = S(2);
S3 = S(3);
S4 = S(4);
S5 = S(5);
S6 = S(6);
% stores at previous timestep
t = obj.t; % this time step
if t == 1
S3old=obj.S0(3);
else
S3old = obj.stores(t-1,3);
end
% climate input
climate_in = obj.input_climate(t,:); % climate at this step
P = climate_in(1);
Ep = climate_in(2);
T = climate_in(3);
% fluxes functions
flux_ei = evap_1(S1,Ep,delta_t);
flux_pn = interception_1(P,S1,rho);
flux_ps = snowfall_1(flux_pn,T,ts);
flux_pr = rainfall_1(flux_pn,T,ts);
flux_ms = melt_1(as,tm,T,S2,delta_t);
flux_fs = refreeze_1(af,as,tm,T,S3,delta_t);
flux_gs = melt_2(gmax,S2,delta_t);
flux_mw = saturation_1(flux_pr+flux_ms,S3,the*S2);
flux_ew = excess_1(S3old,the*S2,delta_t);
flux_eq = flux_mw + flux_gs + flux_ew;
flux_fi = infiltration_4(flux_eq,phi);
flux_rh = effective_1(flux_eq,flux_fi);
flux_eps= effective_1(Ep,flux_ei);
flux_et = evap_22(sw,sm,S4,flux_eps,delta_t);
flux_rd = saturation_1(flux_fi,S4,smax);
flux_l = recharge_6(ksat,c,S4,delta_t);
flux_ls = recharge_7(lmax,flux_l);
flux_lf = effective_1(flux_l,flux_ls);
flux_rf = baseflow_1(kf,S5);
flux_rs = baseflow_1(ks,S6);
% stores ODEs
dS1 = P - flux_ei - flux_pn;
dS2 = flux_ps + flux_fs - flux_ms - flux_gs;
dS3 = flux_pr + flux_ms - flux_fs - flux_mw - flux_ew;
dS4 = flux_fi - flux_et - flux_rd - flux_l;
dS5 = flux_lf - flux_rf;
dS6 = flux_ls - flux_rs;
% outputs
dS = [dS1 dS2 dS3 dS4 dS5 dS6];
fluxes = [flux_ei flux_pn flux_ps flux_pr flux_ms...
flux_fs flux_gs flux_mw flux_ew flux_eq...
flux_rh flux_eps flux_et flux_fi flux_rd...
flux_l flux_lf flux_ls flux_rf flux_rs];
end
% STEP runs at the end of every timestep.
function [out] = evap_1(S,Ep,dt)
%evap_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Evaporation at the potential rate
% Constraints: f <= S/dt
% @(Inputs): S - current storage [mm]
% Ep - potential evaporation rate [mm/d]
% dt - time step size
out = min(S/dt,Ep);
end
function [out] = interception_1(In,S,Smax,varargin)
%interception_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Interception excess when maximum capacity is reached
% Constraints: -
% @(Inputs): In - incoming flux [mm/d]
% S - current storage [mm]
% Smax - maximum storage [mm]
% varargin(1) - smoothing variable r (default 0.01)
% varargin(2) - smoothing variable e (default 5.00)
if size(varargin,2) == 0
out = In.*(1-smoothThreshold_storage_logistic(S,Smax));
elseif size(varargin,2) == 1
out = In.*(1-smoothThreshold_storage_logistic(S,Smax,varargin(1)));
elseif size(varargin,2) == 2
out = In.*(1-smoothThreshold_storage_logistic(S,Smax,varargin(1),varargin(2)));
end
end
function [out] = snowfall_1(In,T,p1,varargin)
%snowfall_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Snowfall based on temperature threshold
% Constraints: -
% @(Inputs): p1 - temperature threshold below which snowfall occurs [oC]
% T - current temperature [oC]
% In - incoming precipitation flux [mm/d]
% varargin(1) - smoothing variable r (default 0.01)
if size(varargin,2) == 0
out = In.*(smoothThreshold_temperature_logistic(T,p1));
elseif size(varargin,2) == 1
out = In.*(smoothThreshold_temperature_logistic(T,p1,varargin(1)));
end
end
function [out] = rainfall_1(In,T,p1,varargin)
%rainfall_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Rainfall based on temperature threshold
% Constraints: -
% @(Inputs): p1 - temperature threshold above which rainfall occurs [oC]
% T - current temperature [oC]
% In - incoming precipitation flux [mm/d]
% varargin(1) - smoothing variable r (default 0.01)
if size(varargin,2) == 0
out = In.*(1-smoothThreshold_temperature_logistic(T,p1));
elseif size(varargin,2) == 1
out = In.*(1-smoothThreshold_temperature_logistic(T,p1,varargin(1)));
end
end
function [out] = melt_1(p1,p2,T,S,dt)
%melt_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Snowmelt from degree-day-factor
% Constraints: f <= S/dt
% @(Inputs): p1 - degree-day factor [mm/oC/d]
% p2 - temperature threshold for snowmelt [oC]
% T - current temperature [oC]
% S - current storage [mm]
% dt - time step size [d]
out = max(min(p1*(T-p2),S/dt),0);
end
function [out] = refreeze_1(p1,p2,p3,T,S,dt)
%refreeze_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Refreezing of stored melted snow
% Constraints: f <= S/dt
% @(Inputs): p1 - reduction fraction of degree-day-factor [-]
% p2 - degree-day-factor [mm/oC/d]
% p3 - temperature threshold for snowmelt [oC]
% T - current temperature [oC]
% S - current storage [mm]
% dt - time step size [d]
out = max(min(p1*p2*(p3-T), S/dt), 0);
end
function [out] = melt_2(p1,S,dt)
%melt_2
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Snowmelt at a constant rate
% Constraints: f <= S/dt
% @(Inputs): p1 - melt rate [mm/d]
% S - current storage [mm]
% dt - time step size [d]
out = min(p1,S/dt);
end
function [out] = saturation_1(In,S,Smax,varargin)
%saturation_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Saturation excess from a store that has reached maximum capacity
% Constraints: -
% @(Inputs): In - incoming flux [mm/d]
% S - current storage [mm]
% Smax - maximum storage [mm]
% varargin(1) - smoothing variable r (default 0.01)
% varargin(2) - smoothing variable e (default 5.00)
if size(varargin,2) == 0
out = In.*(1-smoothThreshold_storage_logistic(S,Smax));
elseif size(varargin,2) == 1
out = In.*(1-smoothThreshold_storage_logistic(S,Smax,varargin(1)));
elseif size(varargin,2) == 2
out = In.*(1-smoothThreshold_storage_logistic(S,Smax,varargin(1),varargin(2)));
end
end
function [out] = excess_1(So,Smax,dt)
%excess_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Storage excess when store size changes (returns flux [mm/d])
% Constraints: f >= 0
% @(Inputs): So - 'old' storage [mm]
% Smax - 'new' maximum storage [mm]
% dt - time step size [d]
out = max((So-Smax)/dt,0);
end
function [out] = infiltration_4(fluxIn,p1)
%infiltration_4
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Constant infiltration rate
% Constraints: f <= fin
% @(Inputs): p1 - Infiltration rate [mm/d]
% fin - incoming flux [mm/d]
out = min(fluxIn,p1);
end
function [out] = effective_1(In1,In2)
%effective_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: General effective flow (returns flux [mm/d])
% Constraints: In1 > In2
% @(Inputs): In1 - first flux [mm/d]
% In2 - second flux [mm/d]
out = max(In1-In2,0);
end
function [out] = evap_22(p1,p2,S,Ep,dt)
%evap_22 3-part piece-wise evaporation
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Threshold-based evaporation rate
% Constraints: f <= S/dt
% @(Inputs): p1 - wilting point [mm]
% p2 - 2nd (lower) threshold [mm]
% S - current storage [mm]
% Ep - potential evapotranspiration rate [mm/d]
% dt - time step size [d]
out = min(S/dt,max(0,min((S-p1)./(p2-p1).*Ep,Ep)));
end
function [out] = recharge_6(p1,p2,S,dt)
%recharge_6
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Recharge to fulfil evaporation demand if the receiving
% store is below a threshold
% Constraints: f <= S/dt
% S >= 0 prevents complex numbers
% @(Inputs): p1 - time coefficient [d-1]
% p2 - non-linear scaling [mm]
% S - current storage [mm]
% dt - time step size [d]
out = min(max(S/dt,0),p1.*max(S,0).^p2);
end
function [out] = recharge_7(p1,fin)
%recharge_7
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Constant recharge limited by incoming flux
% Constraints: -
% @(Inputs): p1 - maximum recharge rate [mm/d]
% fin - incoming flux [mm/d]
out = min(p1,fin);
end
function [out] = baseflow_1(p1,S)
% baseflow_1
% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Flux function
% ------------------
% Description: Outflow from a linear reservoir
% Constraints: -
% @(Inputs): p1 - time scale parameter [d-1]
% S - current storage [mm]
out = p1.*S;
end
function [out] = smoothThreshold_storage_logistic(S,Smax,r,e)
%smoothThreshold_storage_logistic Logisitic smoother for storage threshold functions.
% Copyright (C) 2018 Wouter J.M. Knoben
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Smooths the transition of threshold functions of the form:
%
% Q = { P, if S = Smax
% { 0, if S < Smax
%
% By transforming the equation above to Q = f(P,S,Smax,e,r):
% Q = P * 1/ (1+exp((S-Smax+r*e*Smax)/(r*Smax)))
%
% Inputs:
% S : current storage
% Smax : maximum storage
% r : [optional] smoothing parameter rho, default = 0.01
% e : [optional] smoothing parameter e, default 5
%
% NOTE: this function only outputs the multiplier. This needs to be
% applied to the proper flux utside of this function.
%
% NOTE: can be applied for temperature thresholds as well (i.e. snow
% modules). This simply means that S becomes T, and Smax T0.
% Check for inputs and use defaults if not provided
% NOTE: this is not very elegant, but it is more than a factor 10 faster then:
% if ~exist('r','var'); r = 0.01; end
% if ~exist('e','var'); e = 5.00; end
if nargin == 2
r = 0.01;
e = 5.00;
elseif nargin == 3
r = r{1};
e = 5.00;
elseif nargin == 4
r = r{1};
e = e{1};
end
% Calculate multiplier
Smax = max(Smax,0); % this avoids numerical instabilities when Smax<0
if r*Smax == 0
out = 1 ./ (1+exp((S-Smax+r*e*Smax)/(r)));
else
out = 1 ./ (1+exp((S-Smax+r*e*Smax)/(r*Smax)));
end
end
function [out] = smoothThreshold_temperature_logistic(T,Tt,r)
%smoothThreshold_temperature_logistic Logisitic smoother for temperature threshold functions.
% Copyright (C) 2018 Wouter J.M. Knoben
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.
% Smooths the transition of threshold functions of the form:
%
% Snowfall = { P, if T < Tt
% { 0, if T >= Tt
%
% By transforming the equation above to Sf = f(P,T,Tt,r):
% Sf = P * 1/ (1+exp((T-Tt)/r))
%
% Inputs:
% P : current precipitation
% T : current temperature
% Tt : threshold temperature below which snowfall occurs
% r : [optional] smoothing parameter rho, default = 0.01
%
% NOTE: this function only outputs the multiplier. This needs to be
% applied to the proper flux utside of this function.
% Check for inputs and use defaults if not provided
% NOTE: this is not very elegant, but it is more than a factor 10 faster then:
% if ~exist('r','var'); r = 0.01; end
if nargin == 2
r = 0.01;
end
% Calculate multiplier
out = 1 ./ (1+exp((T-Tt)/(r)));
end
Note
This standardized MARRMoT structure is not identical to the original named model. p01... follow exact source order; s01... are fixed initial stores, defaulting to zero. Octave + optim are required. Solver behavior can differ across runtime versions.
Simulation
from basinforge import Basin, get_model
basin = Basin.from_csv("basin.csv", basin_id="A", area_km2=1200,
q_unit="m3/s", timestep="daily")
q_mm = get_model("MARRMOT_44").simulate(basin)
q_m3s = basin.to_m3s(q_mm)
Supply your actual data and catchment area; temperature-dependent models require a temperature column. For non-daily models, choose an appropriate warmup in model steps.
See calibration, input requirements, sources and verification limitations.