HBV EM 1.0 POPIS METOD & VALIDACE
-1-
OBSAH 1 METODY VÝPOČTU ...................................................................................................... 3 1.1 SRÁŽKY A TEPLOTA VZDUCHU........................................................................... 4 1.2 MODUL SNĚHU.................................................................................................. 5 1.2.1 Klasická metoda teplotního indexu........................................................... 5 1.2.2 Rozšířená metoda teplotního indexu ........................................................ 7 1.2.3 Metoda energetické bilance ................................................................... 10 1.3 PŮDNÍ MODUL ................................................................................................ 17 1.3.1 Potenciální evapotranspirace v klasické metodě TI ................................ 17 1.3.2 Potenciální evapotranspirace v rozšířené metodě TI a energetické bilanci 18 1.4 ODTOK ............................................................................................................. 20 1.5 KALIBRACE – GENETICKÝ ALGORITMUS .......................................................... 21 1.5.1 Reprezentace parametrů ........................................................................ 21 1.5.2 Vytvoření populace ................................................................................. 21 1.5.3 Výběr jedinců a křížení ............................................................................ 22 2 VALIDACE MODELU HBV EM 1.0 ............................................................................... 23 2.1 POROVNÁNÍ MODELŮ DLE OBJEKTIVNÍCH FUNKCÍ ........................................ 25 2.1.1 Černá Nisa ............................................................................................... 25 2.1.2 Kamenice ................................................................................................. 26 2.1.3 Černá Desná ............................................................................................ 26 2.1.4 Vydra ....................................................................................................... 27 2.2 KOMENTÁŘ...................................................................................................... 29 2.2.1 Validace na základě shody měřeného a simulovaného průtoku ............ 29 2.2.2 Problematika simulace vodní hodnoty sněhu......................................... 31 2.2.3 Teplotní index vs. energetická bilance .................................................... 33 3 ZDROJE ...................................................................................................................... 34
-2-
1 METODY VÝPOČTU Model HBV EM 1.0 vychází z konceptuálního modelu HBV ETH, který byl vyvinut na ETH Zürich za účelem modelování srážko-odtokových poměrů v alpských povodích, přičemž byl brán zřetel právě na modul sněhu a ledu (RENNER ET BRAUN, 1990). V roce 1999 byla na Bavorské akademii věd sestavena verze HBV3-ETH9 (KFG, 1999). Model HBV-ETH je založen na řešení rovnice vodní bilance, kterou lze popsat vztahem: 𝑑
𝑃 − 𝐸 − 𝑄 = 𝑑𝑡 (𝑆𝑃 + 𝑆𝑆𝑀 + 𝑆𝑈𝑍 + 𝑆𝐿𝑍) kde je: P srážky [mm] E evapotranspirace [mm] Q odtok [mm] SP zásoba vody ve sněhové pokrývce [mm] SSM zásoba půdní vlhkosti [mm] SUZ zásoba vody v horní nádrži [mm] SLZ zásoba vody ve spodní nádrži [mm]
(1)
Časový krok modelu HBV-ETH je 1 den. Do popisu není vzhledem k cílové aplikaci modelu zahrnuta část týkající se ledovců, kde dochází k redukci rychlosti tání. Schéma modelu HBV-ETH, včetně kalibračních parametrů, je na obrázku 1.
Obrázek 1: Struktura modelu HBV-ETH s tučně vyznačenými kalibrovatelnými parametry (HAGG, 2003).
-3-
HBV EM 1.0 obsahuje tři různé metody stanovení tání (modul sněhu) a stanovení potenciální evapotranspirace (půdní modul). Společná pro všechny je prostorová distribuce jednotlivých složek výpočtu a výpočty probíhající v horní a spodní nádrži.
1.1 SRÁŽKY A TEPLOTA VZDUCHU V modelu jsou teploty a srážky definovány na základě teplotního a srážkového gradientu pro každou buňku zvlášť následovně: 𝑇[𝑥, 𝑦, 𝑐] = 𝑇𝑚 (1 + 𝑇𝐺𝑅𝐴𝐷 ∆𝐻 10−2 ) 𝑃[𝑥, 𝑦, 𝑐] = 𝑃𝑚 (1 + 𝑃𝐺𝑅𝐴𝐷 ∆𝐻 10−4 ) kde je: T[x,y,c] Tm TGRAD ∆H P[x,y,c] Pm PGRAD
(2) (3)
teplota vzduchu v buňce na souřadnici x,y v čase c [°C] měřená teplota v čase c [°C] teplotní gradient [°C/100 m] výškový rozdíl mezi buňkou gridu [x,y] a místem měření srážky v buňce na souřadnici x,y v čase c [mm] měřená srážka v čase c procentuální změna srážky s výškou [%/100m]
Stanovení typu srážek je definováno kritickou teplotou T0. Měřené srážky je možné korigovat parametry RCF pro déšť a SCF pro sníh. Tyto kalibrovatelné parametry korigují vztah mezi nadhodnocením či podhodnocením srážek v místě měření a povodím. Korekci lze popsat: 𝑃𝑙𝑖𝑞 = 𝑅𝐶𝐹 𝑃𝑚 𝑃𝑠𝑛𝑜𝑤 = 𝑆𝐶𝐹 𝑃𝑚 kde je: Pliq RCF Pm Psnow SCF
(4) (5)
korigovaná dešťová srážka [mm] parametr korekce deště [-] měřená srážka [mm] korigovaná sněhová srážka [mm] parametr korekce sněhové srážky [-]
-4-
1.2 MODUL SNĚHU V modulu sněhu lze zvolit jednu ze tří metod výpočtu. Princip, na kterém jsou tyto metody založeny, je popsán v následujících kapitolách. 1.2.1 KLASICKÁ METODA TEPLOTNÍHO INDEXU Modelování tání sněhu je založeno na principu metody teplotního indexu (Degree-Day metoda). Denní úhrny tání jsou tedy odvozovány od teploty a odpovídajícího teplotního indexu, který definuje množství roztátého sněhu na jeden stupeň, pokud teplota překročí definovanou kritickou hodnotu T0. Úhrn tání popisuje rovnice: > 𝑇0 TM{ ≤ 𝑇0 kde je: CMELT CMFi TM T0
𝐶𝑀𝐸𝐿𝑇 = 𝐶𝑀𝐹𝑖 (𝑇𝑀 − 𝑇0) 𝐶𝑀𝐸𝐿𝑇 = 0
(6)
úhrn tání [mm.den-1] teplotní index i-tého dne [mm.°C-1.den-1] denní průměrná teplota vzduchu [°C] kritická teplota, při níž začíná tát sníh [°C]
Hodnota CMELT se liší pro jednotlivé třídy expozice. Pro jižní orientaci je CMELT násobena koeficientem zohledňujícím expozici REXP, pro expozici severní je hodnota CMELT tímto koeficientem dělena. Ostatní expozice zůstávají beze změny pouze jako CMELT (viz. obrázek 2). Pro REXP = 1 model tedy mezi jednotlivými třídami orientace rozdíl tání nezahrnuje. Zohlednění proměnlivé rychlosti tání během je provedeno proměnlivým teplotním indexem (CMF) (viz. obrázek 3).
Obrázek 2: Schématické znázornění vlivu orientace svahu na tání a odtok při zvoleném parametru zohledňujícím expozici REXP = 2 (KONZ, 2003).
-5-
Teplotní index [mm.d-1°C-1]
12
CMAX
10 8 6 4
CMIN
2 0 01/11 01/12 01/1
01/2
01/3
01/4
01/5
01/6
01/7
01/8
01/9 01/10
Datum
Obrázek 3: Příklad průběhu hodnot teplotního indexu během hydrologického roku pro zvolené parametry CMIN = 2, CMAX = 10.
Model zahrnuje tento faktor jako kosinusoidu, přičemž maximální hodnoty nabývá 21. června a minimální 21. prosince. Vývoj CMF v čase během jednoho hydrologického roku lze popsat: 𝐶𝑀𝐹𝑖 = kde je: CMF CMAX CMIN MaxTag i
𝐶𝑀𝐴𝑋−𝐶𝑀𝐼𝑁 2
2𝜋
𝑐𝑜𝑠 (𝑀𝑎𝑥𝑇𝑎𝑔 𝑖) + 𝐶𝑀𝐼𝑁 +
𝐶𝑀𝐴𝑋−𝐶𝑀𝐼𝑁 2
(7)
Degree-Day faktor i-tého dne [mm.°C-1.den-1] maximální Degree-Day faktor [mm.°C-1.den-1] minimální Degree-Day faktor [mm.°C-1.den-1] počet dní v roce minus 1 pořadové číslo dne (i = 0 21. června)
Proces, kdy voda z tajícího sněhu proniká do nižších vrstev sněhové pokrývky, kde opět zamrzá a zůstává, je v modelu rovněž zohledněn. Váhu tohoto procesu definuje parametr CRFR, neboli koeficient „opětovného zamrzání“. Nastává jen v případě, kdy je teplota vzduchu menší než T0. Lze ho popsat jako: 𝐶𝑅𝐸𝐹 = (𝑇𝑀 − 𝑇0) 𝐶𝑀𝐹𝑖 𝐶𝑅𝐹𝑅 kde je: CRFR
(8)
koeficient opětovného zamrzání [-]
Maximální množství vody, které může být ve sněhu tímto způsobem uloženo, je definováno parametrem SLIQMAX následujícím vztahem: 𝑆𝐿𝐼𝑄𝑀𝐴𝑋 = 𝐶𝑊𝐻 𝑆𝑃
(9)
kde je:
-6-
CWH SP
relativní množství vody udržitelné ve sněhové pokrývce [-] vodní hodnota sněhové pokrývky [mm]
Celková bilance vodní hodnoty sněhu SP v čase i je potom určena vztahem: ≤ 𝑇0 𝑆𝑃𝑖 = 𝑆𝑃𝑖−1 + 𝑃𝑆𝑁𝑂𝑊𝑖 𝑇𝑀 { > 𝑇0 𝑆𝑃𝑖 = 𝑆𝑃𝑖−1 − (𝐶𝑀𝐸𝐿𝑇𝑖 − 𝑆𝐿𝐼𝑄𝑖 ) kde je: TM T0 SP PSnow CMELT SLIQ
(10)
teplota vzduchu [°C] kritická teplota vzduchu [°C] vodní hodnota sněhu [mm] vodní hodnota nového sněhu [mm] úhrn tání vodní kapacita sněhu o maximální hodnotě SLIQMAX [mm]
1.2.2 ROZŠÍŘENÁ METODA TEPLOTNÍHO INDEXU Nároky na vstupní data jsou oproti předchozí variantě rozšířeny o pokryv. Parametr REXP, definující rozdíl v tání pro jednotlivé třídy expozice, je nahrazen relativní hodnotou potenciální solární radiace v poměru k horizontální rovině, měnící se v závislosti na čase, zeměpisné šířce, orientaci a sklonu svahů. Potenciální solární radiace, dopadající na horizontální rovinu, je stanovena dle rovnice (Lee, 1966 in DeWalle et Rango, 2008): 𝑅𝑒𝐻 = (𝑆´/𝑒)𝑐𝑜𝑠𝑍 kde je: ReH S´ e Z
(11)
potenciální solární radiace na horizontální rovinu [kJ m-2] solární konstanta = 4896 kJ h-1 m-2 relativní vzdálenost Země/Slunce [-] zenitový úhel [rad]
Relativní vzdálenost Země-Slunce lze určit z rovnice: 2𝜋
𝑒 = (1 + 0,033 𝑐𝑜𝑠 (365 𝐽))−1 kde je: e J
(12)
relativní vzdálenost Země - Slunce [-] pořadové číslo dne v roce (1. leden = 1)
Zenitový úhel Z lze stanovit z rovnice: 𝑍 = 𝑎𝑟𝑐𝑐𝑜𝑠(𝑠𝑖𝑛𝜃 𝑠𝑖𝑛𝛿 + 𝑐𝑜𝑠𝜃 𝑐𝑜𝑠𝛿 𝑐𝑜𝑠𝜔𝑡) kde je:
-7-
(13)
Z θ δ ω t
zenitový úhel [rad] zeměpisná šířka [rad] sluneční deklinace [rad] úhlová rychlost rotace Země = 2π/24 [rad h-1] sluneční čas [h]
Sluneční deklinaci lze určit ze vztahu: 2𝜋
𝛿 = 0,4102 𝑠𝑖𝑛 (365 (𝐽 − 80)) kde je: δ J
(14)
sluneční deklinace [rad] pořadové číslo dne v roce (1. leden = 1)
Pro stanovení denní potenciální solární radiace je nutné určit čas východu a západu Slunce, respektive délku světelného dne. Za předpokladu, že v době východu a západu je zenitový úhel Z = π/2, tedy cosZ = 0 a čas východu a západu je symetricky rozdělen kolem slunečního poledne, lze určit čas ze vztahu: 𝑡1 = −arccos(−𝑡𝑎𝑛𝜃 𝑡𝑎𝑛𝛿) /𝜔 𝑡2 = arccos(−𝑡𝑎𝑛𝜃 𝑡𝑎𝑛𝛿) /𝜔 kde je: t1 t2 θ δ ω
(15) (16)
čas východu Slunce [h] čas západu Slunce [h] zeměpisná šířka [rad] sluneční deklinace [rad] úhlová rychlost rotace Země = 2π/24 [rad h-1]
Pokud je znám čas východu a západu Slunce, lze určit celkovou denní potenciální radiaci na horizontální rovinu podle vztahu: 1
𝑅𝑒𝐻 = (𝑆´/𝑒 2 ) [(𝑡2 − 𝑡1 )𝑠𝑖𝑛𝜃 𝑠𝑖𝑛𝛿 + 𝜔 𝑐𝑜𝑠𝜃 𝑐𝑜𝑠𝛿(𝑠𝑖𝑛𝜔𝑡2 − 𝑠𝑖𝑛𝜔𝑡1 ) kde je: ReH S´ t1 t2 θ δ ω
potenciální solární radiace na horizontální rovinu [kJ m-2] solární konstanta = 4896 kJ h-1 m-2 čas východu Slunce [h] čas západu Slunce [h] zeměpisná šířka [rad] sluneční deklinace [rad] úhlová rychlost rotace Země = 2π/24 [rad h-1]
-8-
(17)
Potenciální solární radiace na svahu o určité orientaci lze stanovit obdobně. Je potřeba nalézt ekvivalentní horizontální povrch, na který dopadá záření ve stejném úhlu jako na vyšetřovaný svah. Výslednou rovnici lze zapsat jako: 1
𝑅𝑒𝑆 = (𝑆´/𝑒 2 ) [(𝑡2 − 𝑡1 )𝑠𝑖𝑛𝜃´ 𝑠𝑖𝑛𝛿 + 𝜔 𝑐𝑜𝑠𝜃´ 𝑐𝑜𝑠𝛿(𝑠𝑖𝑛𝜔𝑡´2 − 𝑠𝑖𝑛𝜔𝑡´1 ) kde je: ReS θ´ ωt´
(18)
potenciální solární radiace na svah [kJ m-2] zeměpisná šířka ekvivalentní roviny [rad] hodinový úhel ekvivalentní roviny [rad]
Zeměpisná šířka a hodinový úhel ekvivalentní roviny lze určit ze vztahů: 𝜃´ = 𝑎𝑟𝑐𝑠𝑖𝑛(𝑠𝑖𝑛𝑘𝑠 𝑐𝑜𝑠ℎ 𝑐𝑜𝑠𝜃 + 𝑐𝑜𝑠𝑘𝑠 𝑠𝑖𝑛𝜃) 𝜔𝑡´ = 𝜔𝑡 + 𝑎 𝑎 = 𝑎𝑟𝑐𝑡𝑎𝑛(𝑠𝑖𝑛ℎ 𝑠𝑖𝑛𝑘𝑠 ) /(𝑐𝑜𝑠𝑘𝑠 𝑐𝑜𝑠𝜃 − 𝑐𝑜𝑠ℎ 𝑠𝑖𝑛𝑘𝑠 𝑠𝑖𝑛𝜃) kde je: ks h a
(19) (20) (21)
sklon terénu [rad] orientace terénu [rad od severu] rozdíl v zeměpisné délce mezi ekvivalentní rovinou a svahem
Východ a západ Slunce jsou stanoveny stejně jako v případě horizontální roviny, do rovnice ale vstupuje místo zeměpisné šířky svahu, zeměpisná šířka ekvivalentní roviny, takže: 𝑡´ = 𝑎𝑟𝑐𝑐𝑜𝑠(−𝑡𝑎𝑛𝜃´ 𝑡𝑎𝑛𝛿) /𝜔 kde je: t´ θ´ δ ω
(22)
čas východu Slunce nad ekvivalentní horizontální rovinou [h] zeměpisná šířka ekvivalentní roviny [rad] sluneční deklinace [rad] úhlová rychlost rotace Země = 2π/24 [rad h-1]
Takto získané hodnoty je nutné porovnat s východem a západem Slunce nad horizontální rovinou o stejné zeměpisné šířce. Výsledné hodnoty t1´ a t2´ pro svah jsou definovány: t1´ = max(t1´, t1) t2´ = min(t2´, t1)
(23) (24)
Výsledná hodnota relativní potenciální solární radiace RR [-] vstupující do modelu místo parametru REXP je určena podílem radiace na rovinu a na svah jako: 𝑅
𝑅𝑅 = 𝑅 𝑒𝑆
(25)
𝑒𝐻
-9-
Pokryv, reprezentován vrstvou Les, dělí území na lesní plochy (hodnota 1) a bezlesí (hodnota 0). Tání je stanoveno opět metodou teplotního indexu, jehož hodnota je v čase proměnná podle stejného principu jako u původního modelu HBV-ETH. Vedle času je však hodnota teplotního indexu závislá právě na vegetačním pokryvu a na výskytu dešťových srážek. Teplotní index je pro každou buňku gridu v čase upravován dle těchto vztahů: a) pokud je PLIQ = 0 1 𝐶𝑀𝐹𝑣 = 𝐶𝑀𝐹 𝑘𝑓 𝐿𝑒𝑠 = { 0 𝐶𝑀𝐹𝑣 = 𝐶𝑀𝐹
(26)
b) pokud je PLIQ > 0 1 𝐶𝑀𝐹𝑣 = 𝐶𝑀𝐹 𝑘𝑓 + 𝐶𝑀𝐹 ∗ 𝑃𝐿𝐼𝑄 ∗ 𝑘𝑓𝑝 𝐿𝑒𝑠 = { 0 𝐶𝑀𝐹𝑣 = 𝐶𝑀𝐹 + 𝐶𝑀𝐹 ∗ 𝑃𝐿𝐼𝑄 ∗ 𝑘𝑐𝑝 kde je: PLIQ Les CMF CMFv kf kfp kcp
(27)
úhrn dešťových srážek [mm] vrstva LandUse reprezentující lesní pokryv/bezlesí hodnota teplotního indexu stanoveného dle rovnice (3.7) [mm °C-1 den-1] hodnota teplotního indexu účastnící se výpočtu [mm °C-1 den-1] koeficient vyjadřující vliv lesního pokryvu na teplotní index [-] koeficient vyjadřující vliv kapalných srážek na teplotní index v lese [-] koeficient vyjadřující vliv kapalných srážek na teplotní index v lese [-]
Koeficienty kf, kfp a kcp mohou být uživatelem v modelu dále upravovány, nikoli však kalibrovány. Doporučené intervaly hodnot, vycházející z výsledků práce Kuusista (1980) jsou ‹0,3;0,7› pro kf, ‹0,02;0,06› pro kfp a ‹0,01;0,03› pro kfcp.
1.2.3 METODA ENERGETICKÉ BILANCE V této variantě je zcela nahrazen původní modul sněhu zjednodušenou metodou energetické bilance, kterou popisuje Walter (2005). Jednotlivé složky rovnice energetické bilance jsou stanoveny na základě minimálních a maximálních denních teplot. Z tohoto důvodu je nárok na vstupní data vyšší. Dalším vstupem je průměrná rychlost větru. V následujícím přehledu jsou uvedeny metody stanovení jednotlivých složek energetické bilance, kterou lze napsat jako: 𝜆𝛥𝑆𝑊𝐸 = 𝑆 + 𝐿𝑎 − 𝐿𝑡 + 𝐻 + 𝐸 + 𝐺 + 𝑃 − 𝑆𝑊𝐸(𝐶𝛥𝑇𝑠 ) kde je: λ ∆SWE S La
skupenské teplo tání ledu = 3,34 x 105 kJ m-3 změna vodní zásoby sněhu [mm] solární radiace [kJ m-2] atmosférická dlouhovlnná radiace [kJ m-2]
- 10 -
(28)
Lt H E G P SWE(C∆Ts)
terestriální dlouhovlnná radiace [kJ m-2] výměna zjevného tepla [kJ m-2] energetický tok vaporizace a kondenzace při povrchu [kJ m-2] vedení tepla půdou do sněhové pokrývky [[kJ m-2] teplo dodané deštěm [kJ m-2] změna zásoby tepla ve sněhové pokrývce
1.2.3.1 Solární radiace Solární radiaci S lze stanovit na základě rovnice: 𝑆 = (1 − 𝐴)𝑇𝑡 𝑅𝑒 kde je: S A Tt Re
(29)
solární radiace [kJ m-2] albedo povrchu [-] atmosférická transmisivita [-] potenciální solární radiace [kJ m-2]
Potenciální solární radiace Re je určena s ohledem na orientaci a sklon svahů dle vztahů uvedených v kap. 1.2.2. Albedo povrchu A je vyjádřeno vztahem (US ARMY CORPS OF ENGINEERS, 1960 IN WALTER 2005): 𝐴 −0,35 2,16 0.46
𝑥 𝐴 = 0,35 − (0,35 − 𝐴𝑥 ) 𝑒𝑥𝑝[−0,177 + 𝑙𝑛 ( 𝐴´−0,35 )
kde je: A A´ Ax
]
(30)
albedo povrchu [-] albedo předchozího dne [-] maximální albedo ~ 0.95
Průměrná hodnota albeda po sněžení je definována vztahem (WALTER, 2005): 4𝑅𝑝 𝜌𝑠𝑛
𝐴 = 𝐴𝑥 − (𝐴𝑥 − 𝐴´) 𝑒𝑥𝑝 [− ( kde je: A A´ Rp ρsn
0,12
)]
(31)
albedo povrchu [-] albedo předchozího dne [-] vodní hodnota nového sněhu [m] hustota nového sněhu [kg m-3]
- 11 -
Hustotu nového sněhu lze vyjádřit jako funkci průměrné teploty vzduchu T a dle vztahu (GOODISON ET AL., 1981): 𝜌𝑠𝑛 = 50 + 3,4(𝑇𝑎 + 15) kde je: ρsn Ta
(32)
hustota nového sněhu [kg m-3] průměrná denní teplota vzduchu [°C]
Atmosférickou transmisivitu lze určit ze vztahu (NDLOVU, 1994 IN WALTER, 2005): 𝑇𝑡 = 0.75[1 − 𝑒𝑥𝑝 (− 𝑅
𝐵
𝑒30
kde je: Tt B Re30 Tx Tn
(𝑇𝑥 − 𝑇𝑛 )2 )]
(33)
atmosférická transmisivita [-] empirický koeficient [-] potenciální solární radiace 30 dní před dnem simulace [MJ m-2 d-1] maximální denní teplota vzduchu [°C] minimální denní teplota vzduchu [°C]
Empirický koeficient B je sezónně a geograficky proměnlivý. Generalizací byly vytvořeny dvě varianty stanovení B, a to pro letní období (90 dní před a po letním slunovratu) a pro zimní období (zbytek roku). Varianty koeficientu B jsou definovány vztahem: 𝐵𝑙é𝑡𝑜 = 0,282𝜙 −0,431 𝐵𝑧𝑖𝑚𝑎 = 0,170𝜙 −0,979 kde je: ϕ
(34) (35)
zeměpisná šířka [rad]
Typ vegetačního pokryvu se na celkové energetické bilanci projevuje skrze redukci solární radiace tak, že je v případě lesa hodnota solární radiace rovna nule. 1.2.3.2 Dlouhovlnná radiace Dlouhovlnná radiace je stanovena na základě Stefan-Boltzmannovy rovnice: 𝐿 = 𝜀𝜎𝑇𝐾4 kde je: L ε σ Tk
(36)
dlouhovlnná radiace [kJ m-2] emisivita povrchu [-] Stefan-Boltzmannova konstanta = 4,89 x 10-11 kJ m-2 K-4 den-1 teplota vyzařujícího tělesa [°K]
- 12 -
Pro atmosférickou dlouhovlnnou radiaci La je emisivita povrchu stanovena na základě rovnice (CAMPBELL ET NORMAN, 1998): 𝜀𝑎 = (0,72 + 0,005𝑇𝑎 )(1 − 0,84𝐶) + 0,84𝐶 kde je: Ta C
(37)
průměrná denní teplota vzduchu [°C] poměr oblačnosti
Pro terestriální dlouhovlnnou radiaci Lt nabývá emisivita povrchu hodnoty 0,97 v případě výskytu sněhové pokrývky a 0,95 v období bez výskytu sněhu. 1.2.3.3 Výměna zjevného tepla Výměna zjevného tepla je stanovena jako (WALTER, 2005): 𝑇𝑠 −𝑇𝑎
𝐻 = 𝐶𝑎 (
𝑟ℎ
)
kde je: H Ca Ts Ta rh
(38)
zjevné teplo [kJ m-2] tepelná kapacita vzduchu ~ 1,29 kJ m-3 °C-1 teplota sněhu [°C] průměrná denní teplota vzduchu [°C] resistence tepelné výměny [den m-1]
Resistenci výměny tepla lze určit ze vztahu: 𝑟ℎ =
𝑧 −𝑑+𝑧ℎ 𝑧 −d+𝑧𝑚 ) 𝑙𝑛( 𝑢 )𝑙𝑛( 𝑇 𝑧𝑚
𝑧ℎ
𝑘2𝑢
kde je: rh zu zT d zm k u
1
(39)
84600
resistence výměny tepla [den m-1] výška měření rychlosti větru [m] výška měření teploty vzduchu [m] výška nulové roviny posunutí ~ 0 m parametr drsnosti - pro sníh ~ 0,0002 m von Karmanova Konstanta = 0,41 rychlost větru [m s-1]
1.2.3.4 Tepelný tok související s výparem a kondenzací Denní tepelný tok způsobený vaporizací a kondenzací lze určit ze vztahu: 𝐸 = 𝜆𝑣 (
𝜌𝑠 −𝜌𝑎 𝑟𝑣
)
(40)
- 13 -
kde je: E λv ρs ρa rv
teplo dodané vaporizací a kondenzací [kJ m2] latentní teplo vaporizace = 2500 kJ kg hustota vodních par při povrchu [kg m-3] hustota vodních par vzduchu [kg m-3] resistence výměny par ~ rh
Zjednodušujícím předpokladem je, že hustota vodních par na povrchu je rovna hustotě nasycených vodních par při teplotě sněhu. Protože je mezi minimální teplotou vzduchu a hustotou vodních par silná korelace, lze přibližně stanovit rovnost mezi hodnotou rosného bodu a minimální denní teplotou. Tedy, že hustota vodních par vzduchu je hustota nasycených vodních par při minimální denní teplotě. Hustota nasycených vodních par při teplotě T [°C] je popsána vztahem (ALLEN, 1990): 16,78𝑇−116,8
𝜌𝑜 = 𝑒𝑥𝑝 ( kde je: ρo R T
𝑇+273,3
1
) [(273,15+𝑇)𝑅]
(41)
hustota nasycených vodních par [kg m-3] termodynamická konstanta = 0,4615 kJ kg-1°K-1 teplota vzduchu [°C]
1.2.3.5 Teplo dodané z půdy Vedení tepla půdy do sněhové pokrývky se jeví jako málo významné. Na základě výzkumu US ARMY CORPS OF ENGINEERS (1960) je vedení tepla půdy do sněhové pokrývky definováno jako konstanta G = 173 kJ m-2 den-1. 1.2.3.6 Teplo dodané deštěm Za předpokladu, že má teplota deště stejnou teplotu jako vzduch, a teplo je dodáváno sněhové pokrývce skrze ochlazování dešťových srážek ve sněhu na hodnotu 0 °C, teplo přijaté sněhovou pokrývkou z deště lze stanovit jako (WALTER, 2005): 𝑃 = 𝐶𝑤 𝑅𝑝 𝑇𝑎 kde je: P Cw Rp Ta
(42)
teplo dodané deštěm [kJ m-2] tepelná kapacita vody = 4,2x109 kJ m-3 °C-1 výška deště [m] teplota vzduchu [°C]
1.2.3.7 Energie uložená ve sněhové pokrývce Během období záporné energetické bilance klesá teplota sněhu v závislosti na teplotě vzduchu, energetické bilanci, hustotě a výšce sněhu. Během období nárůstu přijímané
- 14 -
energie teplota sněhu roste na maximální hodnoty 0 °C. Tepelná kapacita sněhu byla stanovena jako konstantní s hodnotou 2,1 kJ.kg-1°C-1 (WALTER, 2005). 1.2.3.8 Výpočet tání Celková energie je dána vztahem: 𝐸𝑠𝑢𝑚 = 𝑆 + 𝐿𝑎 − 𝐿𝑡 + 𝐻 + 𝐸 + 𝐺 + 𝑃 kde je: Esum S La Lt H E G P
(43)
celková energie [kJ m-2] solární radiace [kJ m-2] atmosférická dlouhovlnná radiace [kJ m-2] terestriální dlouhovlnná radiace [kJ m-2] zjevné teplo [kJ m-2] teplo dodané vaporizací a kondenzací [kJ m2] teplo dodané půdou [kJ m-2] teplo dodané deštěm [kJ m-2]
Úhrn tání je potom určen jako: 𝐶𝑀𝐸𝐿𝑇 =
𝐸𝑠𝑢𝑚 − 𝐶𝑠 𝑆𝑊𝐸𝜌𝑤 (−𝑇𝑠 )
kde je: CMELT Esum Cs SWE ρw -Ts λ
(44)
𝜆𝜌𝑤
úhrn tání [m] suma energií pro tání [kJ m-2] tepelná kapacita sněhu ~ 2,1 kJ kg-1°C-1 vodní hodnota sněhu [m] hustota vody = 1000 kg m-3 teplota sněhu [°C] latentní teplo tání ledu = 3.34 kJ m3
Maximální vodní kapacita sněhu je definována jako: 𝑆𝐿𝐼𝑄𝑀𝐴𝑋 = 0,05 𝑆𝑊𝐸 kde je: SLIQMAX SWE
(45)
Maximální vodní kapacita sněhu [m] vodní hodnota sněhu [m]
Kapalná voda zadržená sněhovou pokrývkou může v závislosti na energetické bilanci opět zamrznout. Množství vody, které může zamrznout, je definováno vztahem (WALTER, 2005): 𝐶𝑅𝐸𝐹 =
𝐸𝑠𝑢𝑚
(46)
𝜆𝜌𝑤
- 15 -
kde je: CREF Esum λ ρw
množství kapalné vody, které může zamrznout [m] suma energií pro tání [kJ m-2] latentní teplo tání ledu = 3.34 kJ m3 hustota vody = 1000 kg m-3
Ekvivalentní srážka je stanovena jako: 𝑅𝑆 = 𝑃𝐿𝐼𝑄 + 𝐶𝑀𝐸𝐿𝑇 − 𝑆𝐿𝐼𝑄 kde je: RS PLIQ CMELT SLIQ
(47)
ekvivalentní srážka [mm] úhrn dešťových srážek [mm] úhrn tání [mm] vodní kapacita sněhu [mm]
- 16 -
1.3 PŮDNÍ MODUL Vstupem do půdního modulu je kapalná srážka PLIQ a úhrn tání CMELT ze všech buněk povodí, v modelu je tento součet označen jako ekvivalentní srážka (RS). Od půdního modulu probíhá výpočet celistvě pro celé povodí. Aktuální zásoba vody v půdním profilu je vyjádřená pomocí hodnoty SSM. Maximální hodnota zásoby vody v modelovém půdním profilu je vyjádřená hodnotou FC. Pokud je hodnota SSM rovna hodnotě FC, je půda plně nasycená vodou. Pokud je půda nasycena nad určitou hranici definovanou hodnotou LP, pak není evapotranspirace limitována zásobou vody v půdním profilu a její aktuální intenzita je rovna potenciální evapotranspiraci, v opačném případě je hodnota aktuální evapotranspirace stanovena dle vztahu: 𝐸𝐴 = 𝐸𝑃 kde je: EA EP SSM LP
𝑆𝑆𝑀
(48)
𝐿𝑃
aktuální evapotranspirace [mm.den-1] potenciální evapotranspirace [mm.den-1] stav nasycení vodní kapacity půdy [mm] maximální hodnota SSM, od které platí EA = EP [mm]
V modelu HBV EM jsou použity dva způsoby stanovení potenciální evapotranspirace, jeden pro klasickou metodu teplotního indexu, druhý je společný pro rozšířenou metodu teplotního indexu a energetickou bilanci.
1.3.1 POTENCIÁLNÍ EVAPOTRANSPIRACE V KLASICKÉ METODĚ TI Potenciální evapotranspirace se v čase mění. V případě modelu s klasickou metodou teplotního indexu má její roční chod tvar sinusoidy, přičemž maximální hodnoty nabývá 1. srpna, minimum (0 mm/den) potom 1. února. Maximální hodnota potenciální evapotranspirace je kontrolována parametrem ETMAX. Příklad vývoje potenciální evapotranspirace, jak ji zahrnuje HBV-ETH je na obrázku 4. Stanovení potenciální evapotranspirace je popsáno rovnicí: 2𝜋
𝐸𝑃 = 0,5 𝐸𝑇𝑀𝐴𝑋 (1 + 𝑠𝑖𝑛 (𝑗 𝑀𝑎𝑥𝑇𝑎𝑔)) kde je: EP ETMAX MaxTag j
(49)
potenciální evapotranspirace [mm/den] maximální hodnota potenciální evapotranspirace [mm/den] počet dní v roce minus 1 pořadí dne od 0 do MaxTag, (j = 0 1. května)
- 17 -
Potenciální evapotranspirace [mm.d-1]
7
ETMAX
6 5 4 3 2 1
01/1 14/1 27/1 09/2 22/2 06/3 19/3 01/4 14/4 27/4 10/5 23/5 05/6 18/6 01/7 14/7 27/7 09/8 22/8 04/9 17/9 30/9 13/10 26/10 08/11 21/11 04/12 17/12 30/12
0
Datum
Obrázek 4: Potenciální evapotranspirace během hydrologického roku pro zvolený parametr ETMAX = 6.
1.3.2 POTENCIÁLNÍ EVAPOTRANSPIRACE V ROZŠÍŘENÉ METODĚ TI A ENERGETICKÉ BILANCI Pro tyto dvě metody je potenciální evapotranspirace stanovena dle metodiky Oudina 2010, čímž je eliminován původní kalibrační parametr ETMAX. Pro teploty vzduchu t ≤ 5 °C je EP = 0, pro teplotu vzduchu t > 5 °C je potenciální evapotranspirace určena vztahem: 𝐸𝑃 =
0,408∗ 𝑅𝑒 (𝑡+5)
kde je: EP Re
(50)
100
potenciální evapotranspirace [mm] potenciální solární radiace [MJ m-2 den-1]
Takto stanovená hodnota potenciální evapotranspirace je následovně upravena podle vegetačního krytu dle výsledků práce Johanssona (1993): 𝐿𝑒𝑠 = { kde je: EP EPL EPC
1 𝐸𝑃𝐿 = 𝐸𝑃 ∗ 1,15 0 𝐸𝑃𝐶 = 𝐸𝑃
(51)
potenciální evapotranspirace [mm] potenciální evapotranspirace v lese [mm] potenciální evapotranspirace na pasece [mm]
Množství vody vstupující do horní nádrže je definováno poměrem půdní vlhkosti a maximální vodní kapacitou půdy a koeficientem BETA. Vztah těchto parametrů je znázorněn na obrázku 5. Je popsáno rovnicí: 𝑆𝑆𝑀
𝐷𝑆𝑈𝑍 = 𝑅𝑆 ( 𝐹𝐶 )𝐵𝐸𝑇𝐴
(52)
- 18 -
kde je: DSUZ RS SSM FC BETA
přítok do horní nádrže [mm] ekvivalentní srážka [mm] zásoba vody v půdě [mm] maximální vodní kapacita půdy [mm] koeficient tvaru [-]
Obrázek 5: Vztah mezi půdní vlhkostí SSM, vodní kapacitou půdy FC, koeficientem BETA a přítokem do horní nádrže DSUZ.
Vyprazdňování půdní nádrže SSM v čase i je v případě, že je hodnota SSM nižší než FC (SSMmax = FC) definováno vztahem: 𝑆𝑆𝑀𝑖 = 𝑆𝑆𝑀𝑖−1 + (𝑅𝑆𝑖 − 𝐷𝑆𝑈𝑍𝑖 ) − 𝐸𝐴𝑖 kde je: DSUZ RS SSM EA
(53)
přítok do horní nádrže [mm] ekvivalentní srážka [mm] zásoba vody v půdě [mm] aktuální evapotranspirace [mm.den-1]
- 19 -
1.4 ODTOK Ve třetím modulu, který je tvořen horní a spodní nádrží, dochází k výpočtu odtoku. Z horní nádrže dochází ke dvěma typům odtoku. K rychlému odtoku Q0 a pomalému odtoku Q1. K rychlému odtoku dochází v případě, že kapacita horní nádrže SUZ překročí mezní hodnotu LUZ. Spodní nádrž zajišťuje pomalý odtok prostřednictvím parametru k2. Maximální množství vody, které se infiltruje z horní nádrže do spodní, která se lineárně vyprazdňuje, kontroluje parametr CPERC. Jednotlivé složky odtoku jsou určeny v čase i ze vztahů (řazeny tak, jak jsou v čase určeny): 𝐶𝑃𝐸𝑅𝑖 = 𝑚𝑖𝑛(𝑆𝑈𝑍𝑖 + 𝐷𝑆𝑈𝑍𝑖 , 𝐶𝑃𝐸𝑅𝐶) 𝑆𝑈𝑍𝑖 = 𝑆𝑈𝑍𝑖−1 + 𝐷𝑆𝑈𝑍𝑖 − 𝐶𝑃𝐸𝑅𝑖 > 0 𝑄0 𝑖 = 𝑘0 (𝑆𝑈𝑍𝑖 − 𝐿𝑈𝑍) 𝑆𝑈𝑍𝑖 − 𝐿𝑈𝑍 { ≤ 0 𝑄0 𝑖 = 0 𝑄1 𝑖 = 𝑘1 𝑆𝑈𝑍𝑖 𝑆𝑈𝑍𝑖 = 𝑆𝑈𝑍𝑖 − 𝑄0 𝑖 − 𝑄1 𝑖 𝐶𝑃𝐸𝑅𝑖 = 𝑚𝑖𝑛(𝑆𝑈𝑍𝑖 , 𝐶𝑃𝐸𝑅𝐶) 𝑆𝐿𝑍𝑖 = 𝑆𝐿𝑍𝑖−1 + 𝐶𝑃𝐸𝑅𝑖 𝑄2 𝑖 = 𝑘2 𝑆𝐿𝑍𝑖 𝑆𝐿𝑍𝑖 = 𝑆𝐿𝑍𝑖 − 𝑄2 𝑖 kde je: Q0 Q1 Q2 SUZ LUZ SLZ CPERC CPER k0 k1 k2
(54) (55) (56) (57) (58) (59) (60) (61) (62)
rychlý odtok z horní nádrže [mm.den-1] pomalý odtok z horní nádrže [mm.den-1] odtok ze spodní nádrže [mm.den-1] stav nasycení horní nádrže [mm] mezní hodnota pro odtok Q0 [mm] stav nasycení spodní nádrže [mm] maximální množství vody, které se může účastnit perkolace [mm] aktuální množství vody které se účastní perkolace [mm] parametr lineární nádrže pro odtok Q0 parametr lineární nádrže pro odtok Q1 parametr lineární nádrže pro odtok Q2
Celkový odtok z povodí je stanoven jako součet všech tří složek odtoku, tedy Q0, Q1 a Q2. Při stanovení celkového odtoku z povodí se vychází z předpokladu, že ekvivalentní srážka opouští malá alpská povodí ve stejný den, ve kterém do systému vstoupila.
- 20 -
1.5 KALIBRACE – GENETICKÝ ALGORITMUS Všechny vytvořené varianty je možné automaticky optimalizovat na základě genetického algoritmu. V původním modelu bylo nutné kalibrovat 19 parametrů (parametry týkající se tání ledu nebyly v případě aplikace v podmínkách České republiky do výpočtu zahrnuty), přičemž lze definovat interval, ve kterém bude hodnota parametru hledána. Algoritmem nalezené hodnoty lze v případě potřeby manuálně upravovat.
1.5.1 REPREZENTACE PARAMETRŮ Genom (jedinec) reprezentuje jednu kalibrační sadu parametrů. Reprezentace kalibračních parametrů je realizována zápisem ve dvojkové soustavě. Každému parametru (genu) je vyhrazeno v genomu 10 míst. Celkově je tedy jeden genom reprezentován číselnou řadou o délce počet parametrů krát 10 míst. Příklad části genomu je znázorněn na obrázku 6.
Obrázek 6: Znázornění prvních 40 číslic genomu s uvedenými kalibračními parametry, které jsou touto částí reprezentovány.
Hodnota genu převedená do desítkové soustavy představuje číslo z intervalu ‹1; 1023›. Tímto číslem je definována pozice v intervalu mezi nastavenými hodnotami hledání kalibračního parametru. Pokud je například interval hledání pro parametr RCF definován jako ‹0,1; 2›, a hodnota genu RCF převedená do desítkové soustavy 333 (viz. obrázek 6), pak je výsledná hodnota parametru RCF stanovena jako: (2 - 0,1)/1023 * 333 = 0,618. Tímto způsobem je zajištěna stejná relativní velikost posunu v intervalu doporučených hodnot jednotlivých kalibračních parametrů.
1.5.2 VYTVOŘENÍ POPULACE První populace je vytvořena za užití generátoru náhodných čísel. Pro každý parametr lze pouze nastavit interval, ve kterém se bude jeho hodnota vyskytovat. Dále je nutné nastavit počet jedinců v populaci. Všechny sady parametrů, reprezentované každým jedincem, jsou postupně použity v modelu a jejich kvalita (fitness) je definována na základě hodnoty Nash-Sutcliffova koeficientu (NASH ET SUTCLIFFE, 1970): 𝑅𝑒𝑓𝑓 = 1 − kde je: Qsim Qobs ̅̅̅̅̅̅ Q obs
2 ∑𝑛 𝑖=1[𝑄𝑠𝑖𝑚 (𝑖)−𝑄𝑜𝑏𝑠 (𝑖)] 𝑛 2 ̅̅̅̅̅̅̅ ∑ [𝑄𝑜𝑏𝑠 (𝑖)−𝑄 𝑜𝑏𝑠 ]
(63)
𝑖=1
simulovaný odtok nebo vodní hodnota sněhu měřený odtok nebo vodní hodnota sněhu průměrný měřený odtok nebo vodní hodnota sněhu
- 21 -
1.5.3 VÝBĚR JEDINCŮ A KŘÍŽENÍ Každému jedinci je na základě fitness přiřazená pravděpodobnost výběru P. Na jejím základě je pak dle principu Roulette Wheel Selection realizován výběr, který lze zapsat jako: random = náhodné číslo, kde 0 ≤ random < 1 suma = 0 pro každého jedince J proveď: suma = suma + P(J) pokud je random < suma tak vyber jedince J Takto vybraní jedinci si v případě, že je splněna pravděpodobnost pro křížení, vymění náhodnou část genetického kódu podle jednobodového křížení. V opačném případě vstupují do nové populace jejich kopie. Při tom může dojít k mutaci, která je realizována výměnou 1 za 0 a naopak. Vedle tohoto procesu je v genetickém algoritmu zohledněn elitismus, takže vždy dva nejlepší jedinci přecházejí rovnou do nové generace. Tento proces je ukončen při dosažení definovaného počtu generačních obměn.
- 22 -
2 VALIDACE MODELU HBV EM 1.0 Model, respektive jeho tři varianty, byly testovány na povodích s významným vlivem sněhové pokrývky. Vybrána byla tři povodí v Jizerských horách, kde probíhá monitoring sněhu a na horním povodí Vydry, kde existuje omezená časová řada vodních hodnot sněhu v oblasti Malé Mokrůvky (obrázek 7).
Obrázek 7: Poloha zvolených experimentálních povodí na mapě průměrných ročních srážkových úhrnů České republiky (upraveno z CHMU).
Kalibrace modelů byla provedena na denních průměrných hodnotách průtoků hydrologických let 2001 - 2005, na letech 2006 - 2010 byla provedena validace. Každý model byl na každém povodí optimalizován desetkrát. Těchto deset sad parametrů bylo dále užito pro validaci. Pro všechny modely jsou shodně nastaveny počáteční intervaly hodnot kalibračních parametrů (tabulka 1), přičemž optimalizace verze metody rozšířeného teplotního indexu se neúčastní parametr REXP, ve verzi energetické bilance potom nefigurují ještě parametry CMIN, CMAX, REXP, CWH a CRFR. Automatické optimalizace žádné z variant se neúčastní parametr T0, který byl nastaven fixně na 0 °C. V grafech jsou jednotlivé modely označeny jako: Dis - metoda teplotního indexu Dis+ - rozšířená metoda teplotního indexu Energy - metoda energetické bilance
- 23 -
Parametry genetického algoritmu jsou pro všechny varianty modelů následující: stanovení fitness dle: Nash-Sutcliffův koeficient (Reff) pravděpodobnost křížení: 0,7 pravděpodobnost mutace: 0,01 velikost populace: 20 počet generací: 200 elitismus - počet míst v další generaci: 2 Tabulka 1:: Zvolené intervaly hledání kalibračních parametrů. Parametr RCF
Spodní mez 0,8
Horní mez 2,5
SCF
0,8
2,5
PGRAD
0
20
TGRAD
-1
-0,3
CMIN
1
4
CMAX
2
10
REXP
1
3
CWH
0,01
0,1
CRFR
0,01
5,0
ETMAX
1
25
LP
1
100
FC
1
100
BETA
0,1
10
LUZ
0
200
CPERC
0
50
k0
0
1
k1
0
1
k2
0
1
Pro každou simulaci je vedle Nash-Sutcliffova koeficientu stanovena střední absolutní chyba (Mean Absolute Error - MAE) : 𝑀𝐴𝐸 =
1 𝑛
∑𝑛𝑖=1|𝑄𝑜𝑏𝑠 (𝑖) − 𝑄𝑠𝑖𝑚 (𝑖)|
(64)
kde je: Qobs měřený průtok i-tého pořadí Qsim simulovaný průtok i-tého pořadí
- 24 -
2.1 POROVNÁNÍ MODELŮ DLE OBJEKTIVNÍCH FUNKCÍ Rozdělení hodnot objektivních funkcí modelů je znázorněno kvartilovými grafy (boxploty) při standardním nastavení (IQR = 1,5). Kvartily jsou tři hodnoty proměnné, které rozdělují neklesající řadu hodnot proměnné na čtyři stejně četné části. První, neboli dolní kvartil, odděluje čtvrtinu statistických jednotek s nejnižší hodnotou statistického znaku od tří čtvrtin jednotek s vyšší, popřípadě stejnou hodnotou znaku. Jemu příslušející kumulovaná relativní četnost je tedy 0.25. Obdobně druhý kvartil, tzv. medián s kumulovanou relativní četností 0.5, dělí statistický soubor na dvě stejně velké části, kde v první jsou statistické jednotky s hodnotou statistického znaku menší než medián, v druhé s hodnotou větší nebo rovnou mediánu. Třetí, neboli horní kvartil odpovídá hodnotě s kumulovanou relativní četností 0.75. Kvartily jsou v grafu znázorněny jako spodní a horní hrana obdélníku (1. a 3. kvartil), medián je znázorněn silnou čarou. Horizontální linie mimo obdélník definují nejnižší a nejvyšší hodnotu, která není odlehlým pozorováním, body nad a pod touto hranicí jsou definované jako pozorované hodnoty menší než dolní kvartil - IQR či větší než horní kvartil + IQR.
2.1.1 ČERNÁ NISA Grafy na obrázku 8 znázorňují rozložení hodnot objektivních funkcí za užití boxplotu. Z rozložení Nash-Sutcliffova koeficientu určeného za celé hydrologické roky jak v případě kalibrace tak validace je patrné ve srovnání s ostatními modely dosažení vyšších hodnot v případě modelu Dis+, kdy je medián všech sad v porovnání s ostatními modely nejvyšší. Zároveň je u modelu Dis+ z grafů patrný vyšší rozptyl hodnot objektivních funkcí, a to v případě Reff i MAE. Statisticky významný rozdíl lze na najít mezi modelem Dis+ a oběma ostatními modely v případě hodnot Nash-Sutcliffova koeficientu při validaci jak pro celé hydrologické roky, tak v případě období tání (grafy 3 a 7 na obrázku 8).
- 25 -
Obrázek 8: Rozložení hodnot objektivních funkcí pro jednotlivé varianty modelů za celé hydrologické roky (grafy 1 - 4) a za období 1.1. - 1.6. (grafy 5 - 8) na povodí Černá Nisa.
2.1.2 KAMENICE Výsledky rozložení hodnot objektivních funkcí jednotlivých modelů za celé hydrologické období nevykazují statisticky významné rozdíly, což je zřejmé i z grafů 1 - 4 na obrázku 9. Rozdíly lze nalézt až při zaměření na validace během epizod tání, kdy je významný rozdíl mezi Dis+ a Energy a to jak v hodnotě Reff (graf 7 na obrázku 9) tak MAE (graf 8 na obrázku 9). Porovnáme-li průměrné hodnoty objektivních funkcí mezi modely Dis a Dis+, jsou hodnoty Reff za období validace 0,70 v případě modelu Dis a 0,73 pro Dis+. Průměrné hodnoty MAE jsou 0,143 pro model Dis a 0,131 pro Dis+. Vzhledem k rozložení hodnot objektivních funkcí však nelze tyto rozdíly mezi variantami Dis a Dis+ považovat za statisticky významné.
Obrázek 9: Rozložení hodnot objektivních funkcí pro jednotlivé varianty modelů za celé hydrologické roky (grafy 1 - 4) a za období 1.1. - 1.6. (grafy 5 - 8) na povodí Kamenice.
2.1.3 ČERNÁ DESNÁ Podle rozložení hodnot objektivních funkcí bylo mezi jednotlivými modely na povodí Černé Desné nalezeno více statisticky významných rozdílů (obrázek 10). Při vyhodnocení za celé období kalibrace je statisticky významný rozdíl mezi hodnotou MAE modelů Dis a Energy (graf 2 na obrázku 10). Další rozdíly jsou potom patrné při zaměření na období tání, kdy při kalibraci je statisticky významný rozdíl mezi hodnotou Reff i MAE mezi modely Dis a Energy. Další statisticky významný rozdíl je ve střední hodnotě MAE mezi modely Dis+ a Energy (grafy 5 a 6).
- 26 -
V období validace je potom významný rozdíl v hodnotě Reff mezi modely Dis a Dis+ (graf 7).
Obrázek 10: Rozložení hodnot objektivních funkcí pro jednotlivé varianty modelů za celé hydrologické roky (grafy 1 - 4) a za období 1.1. - 1.6. (grafy 5 - 8) na povodí Černá Desná.
2.1.4 VYDRA Rozdíl oproti předchozím povodím je zřejmý zejména v rozložení hodnot objektivních funkcí (grafy na obrázku 11). Statisticky významný rozdíl je nalezen v kalibračním období u hodnot Reff mezi modely Dis a Dis+, mezi modely Dis+ a Energy i mezi modely Dis a Energy (graf 1), MAE je potom statisticky odlišná mezi modely Dis+ a Energy (graf 2). Střední hodnoty Reff a MAE během kalibrace za období 1.1 - 1.6 se statisticky významně liší ve středních hodnotách současně u modelů Dis a Energy i u modelů Dis+ a Energy (graf 5). V období validace je statisticky významný rozdíl mezi středními hodnotami R eff modelů Dis a Energy a modelů Dis+ a Energy (graf 2). Během období tání byla statisticky významná odlišnost nalezena u Reff mezi modely Dis a Energy, i mezi modely Dis+ a Energy (graf 7). V případě MAE potom pouze mezi modely Dis+ a Energy (graf 8).
- 27 -
Obrázek 11: Rozložení hodnot objektivních funkcí pro jednotlivé varianty modelů za celé hydrologické roky (grafy 1 - 4) a za období 1.1. - 1.6. (grafy 5 - 8) na povodí Vydra.
- 28 -
2.2 KOMENTÁŘ 2.2.1 VALIDACE NA ZÁKLADĚ SHODY MĚŘENÉHO A SIMULOVANÉHO PRŮTOKU Při hodnocení modelů a sestavování pořadí lze vycházet z více kritérií. Pokud by měly být modely hodnoceny pouze na základě výsledků t-testu, tedy že určitý model by byl považován za lepší pouze v případě, že střední hodnoty sady objektivních funkcí se od hodnot jiného modelu statisticky významně liší, bylo by pořadí modelů, pokud je kritériem Reff průtoku za celé období validace (α = 0,05), následující: 1. Dis+ 2. Dis 3. Energy Přičemž model Dis+ dosáhl statisticky lepšího výsledku než obě další varianty na povodí Černé Nisy a Vydry. Rozdíl v pořadí modelů Dis a Energy je způsoben horší simulací modelu Energy na povodí Vydry. Nebýt tedy povodí Vydry, lze považovat modely Dis a Energy z hlediska Reff za shodné. Pokud se zaměříme v období validace pouze na epizody pravděpodobného tání (1.1. 1.6.), je pořadí modelů následující: 1. Dis+ 2. Dis 3. Energy Dis+ dosáhl statisticky lepších výsledků oproti modelu Dis na povodí Černá Nisa a Černá Desná, lepší než model Energy byl na povodí Černá Nisa, Kamenice a Vydra a je tedy oproti oběma ostatním variantám prokazatelně vhodnější. Mnohem menší rozdíl je opět mezi modely Dis a Energy, kdy stejně jako v případě prvního kritéria pro stanovení pořadí, rozhoduje slabá simulace modelu Energy na povodí Vydra. V ostatních případech nelze mezi modely Dis a Energy považovat rozdíly za statisticky významné. Bude-li kritériem hodnota MAE za celé období validace, statisticky významný rozdíl mezi jednotlivými modely nebyl nalezen na žádném povodí. V případě období tání lze pouze na povodí Kamenice a Vydra prohlásit výsledky modelu Dis+ za lepší ve srovnání s výsledky modelu Energy. Jiné, statisticky významné rozdíly v hodnotách MAE mezi modely na povodích nalezeny nebyly. Pokud bude kritériem hodnota objektivních funkcí mezi měřeným a simulovaným odtokem za celé období validace, lze za nejlepší model považovat variantu Dis+ před verzemi Dis a poslední Energy (tabulka 2). Pokud by však byla posuzována kvalita modelů jen na povodích v Jizerských horách, nelze mezi modely Dis a Energy nalézt významný rozdíl. Pokud by byly modely posouzeny na základě průměrné hodnoty Nash-Sutcliffova koeficientu během epizod tání v letech 2006, 2008 a 2009, kdy rok 2006 reprezentuje tání na konci zimy s výraznou akumulací sněhové pokrývky, rok 2008 nevýraznou zimu s několika epizodami tání během celého jejího průběhu a rok 2009 tání v období bez srážkových událostí, byl by na základě hodnot uvedených v tabulce 3 nejlepší model
- 29 -
Dis+. Rozdíl mezi modelem Dis a Energy je opět ovlivněn slabou simulací na povodí Vydra. Pokud by bylo pořadí sestavenou jen na základě jizerskohorských povodí, je model Energy z pohledu hodnot Reff ve srovnání s modelem Dis lepší. Pří zahrnutí Vydry je situace opačná. Protože neexistuje universální objektivní funkce, které definitivně posoudí kvalitu jednotlivých modelů, nelze užitím jedné či dvou objektivních funkcí zcela přesně posoudit absolutní rozdíl v jejich kvalitě. Všechny varianty modelů poskytovaly srovnatelné výsledky. Na povodí Vydry lze popsat horší simulaci modelem Energy, silnou stránkou tohoto modelu však zůstává menší počet kalibračních parametrů a absence nutnosti kalibrace modulu sněhu. Slabší predikce tohoto modelu na povodí Vydra může být způsobena právě nedostatečně reprezentativními klimatickými daty, které poskytuje meteorologická stanice Churáňov, ale které dostatečně necharakterizují podmínky v oblasti horního povodí Vydry. Tabulka 2: Nejlepší hodnoty objektivních funkcí mezi simulovaným a měřeným průtokem pro všechna povodí i varianty modelu. Černá Nisa Kamenice Černá Desná Vydra Objektivní Varianta fce. modelu Hodnota Pořadí Hodnota Pořadí Hodnota Pořadí Hodnota Pořadí Reff
MAE
Dis
0.520
3
0.731
1-2
0.748
3
0.522
2
Dis+
0.560
1
0.731
1-2
0.780
1
0.550
1
Energy
0.524
2
0.730
3
0.772
2
0.440
3
Dis
0.028
2
0.122
1
0.081
1
1.535
3
Dis+
0.025
1
0.123
2
0.090
2
1.388
1
Energy
0.030
3
0.125
3
0.091
3
1.533
2
Tabulka 3: Průměrné hodnoty Reff všech simulací jednotlivých variant modelů během vybraných epizod tání.
Černá Nisa
Kamenice
Černá Desná
Vydra
rok
Dis
Pořadí
Dis+
Pořadí
Energy
Pořadí
2006
0.84
3
0.88
1
0.85
2
2008
0.71
2
0.74
1
0.68
3
2009
0.36
3
0.54
1
0.39
2
2006
0.85
2
0.87
1
0.81
3
2008
0.66
3
0.73
1
0.69
2
2009
0.67
1
0.64
2
0.52
3
2006
0.68
3
0.84
1
0.73
2
2008
0.52
3
0.59
2
0.6
1
2009
0.62
3
0.68
1
0.66
2
2006
0.55
2
0.62
1
0.36
3
2008
0.38
1
0.35
2
0.15
3
2009
0.57
2
0.62
1
0.29
3
- 30 -
2.2.2 PROBLEMATIKA SIMULACE VODNÍ HODNOTY SNĚHU Při porovnání průměrné hodnoty objektivních funkcí mezi měřenou a simulovanou zásobou sněhu (tabulka 4) lze za nejlepší výsledky považovat simulace modelu Dis, kde přitom není zohledněn lesní pokryv, ani srážkové události během období tání. Tabulka 4: Průměrné hodnoty objektivních funkcí všech simulací vodní hodnoty sněhu jednotlivých variant modelů. Černá Nisa Kamenice Černá Desná Vydra Hodnota Pořadí Hodnota Pořadí Hodnota Pořadí Hodnota Pořadí
Objektivní funkce
Varianta modelu Dis
0.660
1
0.560
1
0.370
2
0.690
1
Reff
Dis+ Energy Dis Dis+ Energy
0.640 0.620 52.600 60.330 58.120
2 3 1 3 2
0.397 0.396 68.080 84.050 84.790
2 3 1 2 3
0.180 0.420 89.880 107.440 88.890
3 1 2 3 1
0.680 -2.170 64.210 63.580 137.970
2 3 2 1 3
MAE
Hodnoty uvedené v tabulce 4 jsou vždy průměrem všech profilů na daném povodí. Pro bližší přehled jsou v tabulce 5 uvedeny hodnoty Nash-Sutcliffova koeficientu pro každý sněhoměrný profil zvlášť. V tabulce je u každého profilu uvedena průměrná hodnota Reff všech deseti simulací a nejlepší hodnota Reff na základě shody mezi měřeným a simulovaným vývojem vodní hodnoty sněhu. Ty simulace, které poskytují na některém z profilů nejlepší simulaci vývoje vodní hodnoty sněhu a jsou současně nejlepšími simulacemi na základě shody měřeného a simulovaného odtoku, jsou zvýrazněny šedým polem. Z přehledu v tabulce 5 je zřejmé, že mezi nejlepší simulace odtoku (podle NashSutcliffova kritéria) nemusí znamenat absolutně nejlepší simulaci vývoje vodní hodnoty sněhu. Analýza hodnot Reff v rámci jednotlivých profilů přináší oproti průměrným hodnotám z tabulky 4 jisté zpřesnění v definici kvality simulace jednotlivých modelů. Společným jevem pro všechna jizerskohorská povodí je rozdíl v úspěchu simulace na profilech v bezlesí a na profilech v lese. Při srovnání průměrných hodnot Reff simulací všech variant modelů, je v případě profilů 1 - 4 nejlepší model Dis, zatímco na profilech 5 a 6 je nejlepší model Energy. V případě absolutně nejlepších simulací je na lesních profilech dosaženo nejlepší shody modelem Dis+, poté modelem Energy. Obdobná je situace na povodí Kamenice, kde je patrný rozdíl v kvalitě simulace modelem Dis na profilech bez lesního pokryvu. Na profilech v lese potom simulace všech modelů srovnatelná. Stejný trend lze pozorovat i na povodí Černé Desné, kde ale celkově nejlepších výsledků dosáhl model Energy. Na povodí Vydry jsou o vývoji sněhu k dispozici velmi omezené údaje, nicméně na profilu 1 je lepší výsledek dosažen modelem Dis, na profilu 2 potom modelem Dis+. Stejně jako v případě simulace průtoku je na tomto povodí patrný odstup v simulaci vodní hodnoty sněhu, i když některé ze simulací modelu Energy poskytují relativně kvalitní výsledky simulace vodní hodnoty sněhu i na tomto povodí.
- 31 -
Tabulka 5: Přehled nejlepších a průměrných hodnot Nash-Sutcliffova koeficientu (Reff) simulací vývoje vodní hodnoty sněhu na jednotlivých profilech zájmových povodí.
Povodí
Černá Nisa
Kamenice
Černá Desná
Vydra
Dis+
Energy
# profilu
Dis Prům. Reff
Max. Reff
Prům. Reff
Max. Reff
Prům. Reff
Max. Reff
1
0.60
0.70
0.53
0.63
0.53
0.66
2
-
-
-
-
-
-
3
0.77
0.88
0.69
0.85
0.65
0.78
4
0.80
0.86
0.70
0.90
0.74
0.84
5
0.72
0.86
0.57
0.93
0.76
0.92
6
0.66
0.73
0.70
0.82
0.69
0.78
1
0.47
0.71
0.25
0.56
0.22
0.61
2
0.47
0.69
0.13
0.47
0.07
0.52
3
0.66
0.78
0.70
0.78
0.63
0.77
4
0.66
0.77
0.66
0.76
0.67
0.83
1
0.21
0.52
-0.11
0.27
0.10
0.35
2
0.62
0.78
0.42
0.64
0.65
0.76
3
0.44
0.73
0.15
0.54
0.37
0.56
4
0.61
0.72
0.52
0.73
0.76
0.80
5
0.40
0.59
0.41
0.62
0.76
0.85
1
0.62
0.95
0.26
0.77
0.07
0.64
2
0.77
0.95
0.83
0.95
-0.51
0.62
Na základě údajů z tabulky 5 lze prohlásit, že modely zohledňující lesní pokryv poskytují obecně horší výsledky v bezlesí, zatímco při srovnání s měřenými daty z lesních lokalit poskytují výsledky lepší. Tento jev však může být způsoben několika faktory, které nelze v modelu ve všech případech zcela přesně postihnout, tj. skutečným stavem ploch, charakterizovaných v modelech jako bezlesí a les a reprezentací takové plochy v modelu s ohledem na rozlišení výpočtu. V modelu Energy se rozdíl mezi lesem a bezlesím projevuje redukcí krátkovlnné sluneční radiace na nulovou hodnotu v případě lesa, v modelu Dis+ je tání v lese redukováno koeficientem kf, přičemž nejlepších výsledků bylo dosaženo s hodnotou 0,7. Toto dělení se však jeví jako nedostatečné, zejména v případě oblastí, kde dochází k zarůstání odlesněných ploch.
- 32 -
2.2.3 TEPLOTNÍ INDEX VS. ENERGETICKÁ BILANCE Model s výpočtem tání sněhu na bázi zjednodušené energetické bilance poskytuje ve srovnání s modely teplotního indexu srovnatelně kvalitní výsledky na povodí v Jizerských horách. V simulaci vodní hodnoty sněhu vykazuje model Energy stejný trend jako model Dis+, kdy by bylo vhodné přesněji definovat jednotlivé plochy. Na povodí Černé Desné je odhad vývoje vodní hodnoty sněhu modelem Energy v porovnání s ostatními modely nejpřesnější. Důvodem může být nejmenší vzdálenost povodí od meteorologické stanice Desná-Souš, která je zdrojem vstupních teplot vzduchu. Citlivost na reprezentativnost vstupních dat se projevuje na povodí Vydry. Vzhledem k omezeným údajům o vývoji sněhu na tomto povodí nelze provést kvalitní analýzu simulace vodní hodnoty sněhu, nicméně na základě vyhodnocení shody měřeného a simulovaného průtoku je model Energy z použitých modelů nejslabší, ačkoliv některé simulace poskytují i na tomto povodí relativně spolehlivé výsledky. Důvodem je pravděpodobně výraznější odlišnost klimatických poměrů mezi povodím Vydry a meteorologickou stanicí Churáňov. Silnou stránkou modelu Energy je redukce kalibračních parametrů. Metodu zjednodušené energetické bilance lze tedy doporučit v oblastech s dostatečně reprezentativními zdroji dat.
- 33 -
3 ZDROJE ALLEN, R.G. (ED.) (1990): Evapotranspiration and Irrigation Water Requirements, American Society of Civil Engineering (ASCE), New York, NY, 332 s. CAMPBELL, G.S., NORMAN, J.M., (1998): An Introduction to Environmental Biophysics, second ed., Springer, New York, 286 s. CHMU (2013a): [online], [cit. 2015-01-15] Mapa průměrných ročních úhrnů srážek, vytvořeno 29.1.2013, dostupné z WWW:
DEWALLE, D.R., RANGO, A. (2008): Principles of Snow Hydrology. Cambridge University Press. New York, NY. GOODISON, B.E., FERGUSON, H.L., MCKAY, G.A. (1981): Measurement and data analysis. In: Gray, D.M., Male, D.H. (Eds.), Handbook of Snow, Pergamon Press, Toronto, ON, s. 191–274. HAGG, W. (2003): Auswirkungen von Gletscherschwund auf die Wasserspende hochalpiner Gebiete, Vergleich Alpen – Zentralasien, Münchener Geographische Abhandlungen in Münchener Universitätsschriften, Fakultät für Geowissenschaften, ISBN 3 925308 57 1 JOHANSSON, B. (1999): The relationship between catchment characteristics and the parameters of a conceptual runoff model--A study in the south of Sweden, Contribution to the Second International Conference on FRIEND, October, 1993, Braunschweig, IAHS Publication No, 221, s, 475-482. KONZ, M. (2003): HBV3-ETH9 User Manual, Kommission für Glaziologie der Bayerischen Akademie der Wissenschaften, München, 29. s. LEE, R. (1963): Evaluation of Solar Beam Irradiation as a Climatic Parameter of Mountain Watesheds, Hydrol. Paper. No. 2, Ft. Collin, CO: Colorado State University NDLOVU, L.S., (1994): Weather Data Generation and Its Use in Estimating Evapotranspiration. PhD Thesis. Department of Biosystems Engineering (Engineering Science). Washington State University. Pullman, WA. 143 s. OUDIN, L., MOULIN, L., BENDJOUDI, H., RIBSTEIN, P. (2010): Estimating potential evapotranspiration without continuous daily data: possible errors and impact on water balance simulation. Hydrological Sciences Journal, 55, 2, s. 209 - 222. U.S. ARMY CORPS OF ENGINEERS (1960): Engineering and Design—Runoff from Snowmelt. Engineering Manual 1110-2-1406 WALTER, M.T., BROOKS, E.S., MCCOOL, D.K., KING, L.G., MOLNAU, M., BOLL, J. (2005): Process-based snowmelt modeling: Does it require more input data than temperature-index modeling? Journal of Hydrology 300 (1-4), s. 65-75.
- 34 -
