3D hydrogeological modelling at a mining development site with complex hydrostratigraphy in northern Finland

(1)

3D hydrogeological modelling at a

mining development site with complex hydrostratigraphy in northern Finland

SUSANNE ÅBERG

DEPARTMENT OF GEOSCIENCES AND GEOGRAPHY A95 / HELSINKI 2021 ACADEMIC DISSERTATION

To be presented, with the permission of the Faculty of Science of the University of Helsinki, for public examination in auditorium A110, Chemicum, Kumpula Campus on 22 October 2021, at 10 morning.

(2)

Author´s address: Susanne Åberg

Department of Geosciences and Geography

P.O.Box 64

00014 University of Helsinki Finland

VXVDQQHDEHUJ#KHOVLQNL¿

Supervised by: Docent Kirsti Korkka-Niemi

Department of Geosciences and Geography University of Helsinki

Professor Emeritus Veli-Pekka Salonen Department of Geosciences and Geography University of Helsinki

Reviewed by: Professor Diana Allen Department of Earth Sciences Simon Fraser University Canada

Accociate Professor Martin Ross

Department of Earth and Environmental Sciences

University of Waterloo

Canada

Opponent: Professor Masaki Hayashi Department of Geoscience University of Calgary Canada

(3)

ISSN 1798-7911

ISBN 978-951-51-6590-9 (pbk.) ISBN 978-951-51-6591-6 (PDF) KWWSHWKHVLVKHOVLQNL¿

Punamusta Joensuu 2021

(4)

67°30’02’’ pohjoista leveyttä

ja 26°40’ itäistä pituutta ilman sekstanttia osapuilleen runo Kitisen rannalla kun Kotisuvannosta jäät

kun sauvoimet ja soukat veneet ja tervan haju kun Kotiaavalla kurki

ja villihanhet yli.

…

Pakolliset kuviot

voimalaitokset kaukana vielä monen vuoden takana etelässä:

satatuhatta tai miljoona litraa vettä sekunnissa nivan kohdalla Poronkuvanivan

jättiläistelejä Kitisen rannalla Suurin Kirjaimin.

Nousviikolla kesä.

Excerpts from Veikko Haakana’s poem Jonain iltana yönä aamulla from 1972

(5)

Åberg S., 2021. 3D hydrogeological modelling at a mining development site with complex hydro- stratigraphy in northern Finland3XQDPXVWD-RHQVXXSDJHVDQG¿JXUHV

Abstract

There is an increase in mining projects in northern regions in the areas of complex Quaterna- ry sediments overlaying weathered bedrock, where groundwater systems are understudied.

This study increase the understanding of the interaction of surface and groundwaters in an eco- logically sensitive northern area where mining and other human activities occur.

Groundwater and surface water interactions in riverbanks and open mire areas are critical environments of groundwater-dependent ecosystems (GDEs). Mires are a prevailing feature in northern boreal and subarctic regions, and their interactions with groundwater and other surface water systems are usually poorly understood. If nutrient-rich fens are located in an area where mining activities are planned to take place, a detailed understanding of the water balance and groundwater discharge and recharge patterns are needed, especially how groundwater is linked to properties of the subsurface sediments. Northern areas glaciated and located in weak glacial erosion areas during the last ice ages, could form a complex sedimentary succession interlayered with low-conductivity till and variable sorted sediments having hydraulic conductivity that is orders of magnitudes higher. The glacial–inter- glacial cycles can enable the formation of highly heterogeneous and scattered sediment units, which are challenging to model.

The study area is located in the western part of a Natura 2000-protected Viiankiaapa mire, ZKLFKOLHVDERYHDKLJKJUDGHDQGVLJQL¿FDQW Ni-Cu-PGE deposit. The hydrology of the reg-

XODWHG5LYHU.LWLQHQQHDUE\D൵HFWVWKHZHVWHUQ part of the Viiankiaapa mire, presumable sup- SRUWLQJWKHKDELWDWVRIJURXQGZDWHULQÀXHQFHG mire species. The construction of hydroelectric power plants and river regulation have changed the hydrology of the study area from the 1970s RQZDUGV7KH¿UVWREMHFWLYHRIWKLVVWXG\ZDV to 1) characterize the groundwater recharge/dis- FKDUJHDQGÀRZSDWWHUQVRIJURXQGZDWHULQ4XD- ternary sediments and shallow bedrock, 2) exam- LQHWKHH൵HFWRIK\GURVWUDWLJUDSKLFFRPSOH[LW\

RQJURXQGZDWHUUHFKDUJHGLVFKDUJHDQGÀRZSDW- WHUQVDQGLQYHVWLJDWHWKHH൵HFWRIULYHUUHJX- lation on the hydrology of the Viiankiaapa mire DQGWKHKDELWDWVRIJURXQGZDWHULQÀXHQFHGPLUH species. The second objective of the study was WRFUHDWHDZRUNÀRZIRU'JURXQGZDWHUÀRZ modelling suitable for baseline studies located in weak glacial erosion areas. The 3D groundwater ÀRZPRGHOOLQJZDVDSSOLHGWRPRGHOWKHJURXQG- water recharge/discharge patterns of the study area. Flood modelling was used to model the SUHUHJXODWLRQÀRRGFRYHUDJHRIWKHVWXG\DUHD DQGWRHYDOXDWHLWVH൵HFWRQWKHSUHVHQWKDELWDWV RIJURXQGZDWHULQÀXHQFHGPLUHVSHFLHV6WDEOH isotopes of water and thermal imaging with an unmanned aerial vehicle were used to determine groundwater discharge locations and to verify the JURXQGZDWHUÀRZPRGHOV

7KHUHVXOWVRIJURXQGZDWHUÀRZPRGHOOLQJ indicated that water from the Viiankiaapa mire ÀRZVWRZDUGVWKH5LYHU.LWLQHQDQGGLVFKDUJ- es locally within the mire area and along the shoreline of the river. Groundwater recharge/

(6)

GLVFKDUJH DQG ÀRZ SDWWHUQV ZHUH D൵HFWHG E\

the high complexity of the hydrostratigraphy of VXU¿FLDOGHSRVLWVZHDWKHUHGIUDFWXUHGEHGURFN and small-scale topographical variation within the mire area.

The results indicated that high variation in hydraulic conductivities dispersed recharge/discharge patterns in hydrostratigraphically detailed models compared to more simple models. The ZHDWKHULQJSUR¿OHRIWKHIUDFWXUHGEHGURFNDQG the variation in hydraulic conductivity were also found to be important in modelling the connections of shallow groundwater in sediments and in the topmost part of the bedrock.

Construction of the Matarakoski and Kelu- NRVNLSRZHUSODQWVLQWKH5LYHU.LWLQHQD൵HFW- ed the hydrology of the study area by reducing VSULQJÀRRGVDQGE\UDLVLQJWKHULYHUVWDJH7KH river stage rise caused a reduction in the hydraulic gradient towards the River Kitinen. The reduction in the hydraulic gradient raised the water table, increased groundwater discharge in the

western part of the mire and decreased groundwater discharge into the River Kitinen. Further- more, the modelling results indicated that half of the present habitats of the studied groundwater- LQÀXHQFHGSODQWVSHFLHVRFFXULQDUHDVD൵HFWHG by the regulation of the river.

7KHPDMRUSKDVHVRIWKHFUHDWHGZRUNÀRZ ZHUHGH¿QLWLRQRIWKHK\GURVWUDWLJUDSK\EDVHG RQWKHGL൵HUHQWLDWLRQRIJODFLDOWLOOVDQGLQWHUOD\- ered sorted sediments, 2) the use of models such as MODFLOW-NWT, allowing the modelling of XQFRQ¿QHGVHWWLQJVLQWKHPRGHODQGPRGL-

¿FDWLRQRIWKHK\GURVWUDWLJUDSK\RIWKHPRGHO EDVHGRQFDOLEUDWLRQHYDOXDWLRQDQGYHUL¿FDWLRQ results in an iterative manner. According to the results of this study, constructing a 3D geological PRGHODQGD'JURXQGZDWHUÀRZPRGHOLVJLY- ing valuable support for baseline studies and in the early planning stages of mining. The under- VWDQGLQJRISDVWDQGSUHVHQWDQWKURSRJHQLFLQÀX- ences can be valuable for these baseline studies.

(7)

Tiivistelmä (in Finnish)

Kaivosprojektien määrä on lisääntymässä poh- joisilla alueilla, joissa kompleksiset kvartäärise- dimenttikerrostumat peittävät rapakalliota ja joissa alueen pohjavesivarastot tunnetaan usein hei- kosti. Tämän tutkimuksen tarkoituksena on lisätä ymmärrystä pohjavesien ja pintavesien liikkeistä ja vuorovaikutuksista ekologisesti herkillä poh- joisilla alueilla, joissa on kaivostoimintaa tai sen suunnittelua, taikka muuta ihmisen toimintaa.

Pinta- ja pohjaveden vuorovaikutus on tärkeää tunnistaa, sillä jokien ranta-alueilla ja aa- pasoilla on usein tärkeitä pohjavesivaikutteisten ekosysteemien esiintymisalueita. Aapasuot ovat yleisiä pohjoisella boreaalisella vyöhykkeellä, mutta niiden yhteys muihin vesistöihin ja poh- javesivarastoihin on usein huonosti tunnettu.

Mikäli ravinnerikkaille lettoalueille suunnitel- laan kaivostoimintaa, tarvitaan yksityiskohtais- ta tietoa ja ymmärrystä paikallisesta vesitasees- ta sekä pohjaveden muodostumis- ja purkautu- miskuvioista, jotka linkittyvät suoraan maaperän rakenteeseen. Pohjoiset alueet ovat jäätiköityneet toistuvasti viimeisten jääkausien aikana synnyt- täen monimutkaisen sedimenttisukkesion. Hei- kon jäätikön kulutuksen alueella vuorottelevat al- haisen vedenjohtavuuden moreenit ja vaihtelevat lajittuneet sedimentit, joilla on usein moreeneja selvästi korkeampi hydraulinen johtavuus. Gla- siaali-interglasiaalisyklien synnyttämät vaihtelevat ja hajanaiset sedimenttiyksiköt ovat hyd- rogeologisen mallinnuksen kannalta haastavia ympäristöjä.

Tutkimusalue sijaitsee Natura 2000 -suojel- lun Viiankiaapa nimisen aapasuon alueella, jon- ka alapuolella on rikas Ni-Cu-PGE esiintymä.

Aapasuon vieressä virtaava säännöstelty Kitisen joki vaikuttaa hydrologialtaan Viiankiaavan län- tiseen osaan, jossa esiintyy oletetusti pohjavesivaikutteisia ekosysteemejä. Vesivoimalaitosten

rakentaminen ja jokien säännöstely ovat muuttaneet tutkimusalueen hydrologiaa 1970-luvul- ta lähtien. Tutkimuksen ensimmäinen tavoite oli 1) karakterisoida pohjaveden muodostuminen/

purkautuminen ja virtauskuviot maaperän pohja- vesisysteemissä ja kallion ylimmässä osassa, 2) WXWNLD PDD MD NDOOLRSHUlQ K\GURVWUDWLJUD¿VHQ mallin kompleksisuuden lisäämisen vaikutusta mallinnettuihin pohjaveden muodostumis/purkautumis- ja virtauskuvioihin, ja sekä 3) selvit- tää, miten joen säännöstely on vaikuttanut Viian- kiaavan hydrologiaan ja pohjavesivaikutteisten suokasvilajien esiintymiin. Toinen päätavoite oli kehittää heikon jäätikön kulutuksen alueelle si- joittuville kaivostoiminnan suunnittelualueille soveltuva 3D-pohjaveden virtausmallinnuksen työnkulku. 3D-pohjaveden virtausmallinnus- ta sovellettiin pohjaveden muodostumis/purkautumis- ja virtauskuvioiden mallintamiseen tutkimusalueella. Tulvamallinnusta käytettiin puolestaan säännöstelyä edeltävän tulvan vaiku- tusalueen laajuuden selvittämiseen ja sen arvioi- miseen, miten tulviminen on vaikuttanut suoalu- een pohjavesivaikutteisten kasvien esiintymiin.

Veden hapen ja vedyn stabiileja isotooppeja ja lämpökamerakuvausta miehittämättömällä len- toalustalla hyödynnettiin pohjaveden purkautu- misalueiden tunnistamiseen sekä pohjaveden vir- WDXVPDOOLHQYHUL¿RLPLVHHQ

Pohjaveden virtausmallinnusten tulosten pe- rusteella pohjavesi virtaa tutkimusalueella pääa- siassa suolta kohti Kitistä. Viiankiaavalla osa suovedestä imeytyy pohjavesivyöhykkeeseen ja virtaa kohti Kitistä purkautuen paikoin suolla ja joen rannan lähteissä. Pohjaveden muodostumis/purkautumis- ja virtauskuvioihin vaikuttaa PDDSHUlQNHUURVWHQK\GURVWUDWLJUD¿QHQUDNHQQH kallioperän pintaosan rikkonaisuus ja rapauma VHNlWRSRJUD¿DQYlKlLVHWYDLKWHOXWVXROOD

(8)

Pohjaveden virtausmallien tulosten mukaan maa- ja kallioperän hydraulisen johtavuuden suuri vaihtelu on keskeinen syys siihen, että pohjaveden muodostumis/purkautumiskuvio on yksityiskohtaisempi ja hajaututempi hydrostratigrafisesti yksityiskohtaisemmissa malleissa kuin yksinkertaisemmissa malleissa.

Rapakallion esiintyminen ja luonne vaikuttavat maaperän pohjaveden ja kallion pohjaveden yhteyksiin, sillä rapautuneen kallion hydraulinen johtavuus vaihtelee paljon riippuen sen ra- pautumisasteesta.

Matarakosken ja Kelukosken vesivoimalaitosten rakentaminen on muuttanut tutkimusalueen hydrologiaa vähentäen joen tul- vimista ja nostaen joen vedenpinnan tasoa.

Joen vedenpinnan tason nousu on loiventanut pohjaveden hydraulista gradienttia kohti jokea.

Hydraulisen gradientin muutos on puolestaan nostanut pohjaveden pintaa ja lisännyt pohjaveden purkautumista suon länsilaidalla samalla vähentäen pohjaveden purkautumista Kitiseen.

Mallinnustulokset osoittavat lisäksi, että joen säännöstelyn aiheuttamat muutokset ovat

muuttaneet pohjaveden pintaa alueilla, joissa tutkittuja pohjavesivaikutteisia lajeja esiintyy yleisesti.

Kehitetyn 3D-pohjaveden virtausmallinnuksen työnkulun pääkohtina ovat 1) määrittää PDDSHUlQMDNDOOLRSHUlQK\GURVWUDWLJUD¿VHW\N- siköt perustuen jäätikön synnyttämiin moreeni- yksiköihin ja niiden välissä esiintyviin lajittu- neisiin sedimentteihin sekä kallioperän rikkonai- suusvyöhykkeisiin, 2) käyttää vapaan akviferin mallinnukseen sopivia mallinnuskoodeja kuten MODFLOW-NWT:tä ja 3) muuttaa hydrostrati- JUD¿VWD UDNHQQHWWD SHUXVWXHQ SRKMDYHGHQ YLU- WDXNVHQPDOOLQQXVNDOLEURLQWLMDYHUL¿RLQWLWX- loksiin iteratiivisella tavalla. Tämän tutkimuksen tulokset osoittavat, että 3D-geologinen mal- linnus ja 3D-pohjaveden virtausmallinnus ovat yhdessä käyttökelpoinen työkalu, jota voidaan hyödyntää jo kaivostoiminnan suunnitellun alku- vaiheessa. Pohjaveden ja vesistöjen yhteyksien karakteisoinnilla sekä tutkimalla historiallisen ja nykyisen ihmistoiminnan vaikutusta alueen hydrologiaan ja hydrogeologiaan voidaan saavuttaa arvokasta tietoa jo perustilatutkimusvaiheessa.

(9)

Acknowledgements

I would like to show my greatest appreciation to my supervisors, Associate Research Professor Kirsti Korkka-Niemi, for her generous support and encouragement during this doctoral study, and to Prof. Emer. Veli-Pekka Salonen, who orig- inally suggested this project to me. Without their guidance and persistent help, this thesis would not have been possible. I would like to express my gratitude to all my co-authors: Associate Re- search Professor Kirsti Korkka-Niemi, Dr An- nika Åberg, Dr Anne Rautio, and Prof. Emer.

Veli-Pekka Salonen for their insightful contributions to the manuscripts and helping me to un- derstand the hydrogeological and other aspects of the study area.

Pre-examiners Prof. Diana Allen and Ass.

Prof. Martin Ross are gratefully acknowledged for their constructive criticism and valuable comments on the thesis.

I would like to thank AA Sakatti Mining Oy for the possibility to use their data to enable the construction of the model and interpretation of the results. I would like to thank Dr Anne Rau- tio, Prof. Pertti Lamberg, Jukka Jokela, Laura Paakkonen, Ulla Syrjälä, Janne Siikaluoma and Joanna Kuntonen-van’t Riet for their valuable comments and help with the data.

I am indebted to Rudi-Matti Suoknuuti, Ta- tu Lahtinen, Jyri Laakso, Enni Suonperä, Dr Emilia Koivisto, Mimmi Takalo and Harri Tur- WLDLQHQIRUWKHLUVXSSRUWGXULQJWKH¿HOGZRUN I am deeply grateful to Dr Seija Kultti, Dr Anu Kaakinen, Prof. Emer. Gustaf Olsson, Ass. Prof.

Pietari Skyttä and Ass. Prof. David Whipp for their insightful comments and suggestions for some manuscripts. I would also like to extend my thanks to Roy Siddall for language revision

of all the manuscripts related to this thesis.

,KDYHJUHDWO\EHQH¿WHGIURPWKHGDWDDQG knowledge acquired from Mikko Huokuna, Kai- sa Puolamaa and Taina Kojola from SYKE, Juha- Petri Kämäräinen, Aulis Ruotsalainen and Jari Uusikivi from the Centre for Economic Devel- opment, Transport and the Environment, Juho Päiväniemi from Kemijoki Oy, and local resi- dent Heino Riipi. I would also like to express my gratitude to Richard Winston for guidance with ModelMuse and John Doherty for guidance with PEST.

I am grateful for the support and encouragement of colleagues at the Department of Geosci- ences and Geography, creating a nice and peace- ful working environment. I would especially like to thank Salla Eeva, Tuomas Junna, Juha Järvinen, Dr Yurui Zhang and Peter Howett, who VKDUHGWKHR൶FHZLWKPHGXULQJWKHVH\HDUV, would like to thank my former University of Hel- sinki colleague Dr Juulia-Gabrielle Moreau for friendship and cheering me up during busy times.

I would also like to express my gratitude to Maa- ja Vesitekniikan tuki ry, the K.H. Ren- lund Foundation and the University of Helsinki

*HRGRFIRUWKH¿QDQFLDOVXSSRUWWKDWPDGHWKLV study possible. Additional data collection was performed in several projects at the University, funded by AA Sakatti Mining Oy.

I owe an important debt to many people who helped me and my sister during our adventures LQ&RORUDGRGXULQJDÀRZPRGHOOLQJFRXUVH, would like to thank my parents, my sisters and UHODWLYHVZKRHQFRXUDJHGPHRQP\VFLHQWL¿F MRXUQH\DQGR൵HUHGJRRGGLQQHUV6SHFLDOWKDQNV go to my sister Annika for sharing inspiring ideas and the workload during this study.

(10)

(11)

List of original publications

This thesis is based on the following publications:

I Åberg, S.C., Korkka-Niemi, K., Rautio, A., Salonen, V.-P. & Åberg, A.K. 2019.

Groundwater recharge/discharge patterns and groundwater–surface water interactions in a sedimentary aquifer along the River Kitinen in Sodankylä, northern Fin- land. Boreal Environment Research 24, 155–187.

II Åberg, S.C., Åberg, A.K., & Korkka-Niemi, K. 2021. Three-dimensional hydro- VWUDWLJUDSK\DQGJURXQGZDWHUÀRZPRGHOVLQFRPSOH[4XDWHUQDU\GHSRVLWVDQG weathered/fractured bedrock: evaluating increasing model complexity. Hydroge- ology Journal 29, 1043–1074.

,,, cEHUJ6&.RUNND1LHPL.5DXWLR$ cEHUJ$.7KHH൵HFWRIULYHUUHJ- XODWLRQRQJURXQGZDWHUÀRZSDWWHUQVDQGWKHK\GURORJLFDOFRQGLWLRQVRIDQDDSD mire in northern Finland (Submitted to Journal of Hydrology: Regional Studies).

The publications are referred to in the text by their roman numerals.

Author’s contribution

I The study was designed by K. Korkka-Niemi and S.C. Åberg. S.C. Åberg con- VWUXFWHGWKHJURXQGZDWHUÀRZPRGHOVDQGSUHSDUHGPRVWRIWKH¿JXUHVDQGWKH WDEOHV$5DXWLRSUHSDUHGWKHLVRWRSHUHODWHG¿JXUHVDQG$.cEHUJWKHJHRORJL- FDOPRGHOUHODWHG¿JXUHV6&cEHUJSHUIRUPHGWKHVWDWLVWLFDODQDO\VHVUHODWHGWR groundwater monitoring data. The manuscript was jointly written by S.C. Åberg and K. Korkka-Niemi, with contributions and comments from the co-authors.

II The study was designed by all the authors. S.C. Åberg was responsible for the con- VWUXFWLRQFDOLEUDWLRQDQGYHUL¿FDWLRQRIWKHJURXQGZDWHUÀRZPRGHOV6&cEHUJ SURGXFHGJURXQGZDWHUÀRZPRGHOOLQJUHODWHG¿JXUHVDQGWDEOHV$.cEHUJKDG the main responsibility for the geological/hydrostratigraphic models. However, S.C. Åberg participated in every model designing to a varying degree and was GHVLJQLQJWKHK\GURVWUDWLJUDSKLFXQLWVXVHGIRUWKHJURXQGZDWHUÀRZPRGHOV7KH manuscript was jointly written by all the authors.

III The study was designed by S.C. Åberg and K. Korkka-Niemi. S.C. Åberg con- VWUXFWHGWKHJURXQGZDWHUÀRZPRGHOVDQGÀRRGPRGHOVDQGSUHSDUHGPRVWRIWKH

¿JXUHVDQGWDEOHV$5DXWLRSUHSDUHGWKHLVRWRSHUHODWHG¿JXUHVDQG$.cEHUJ

(13)

WKHJHRORJLFDOPRGHOUHODWHG¿JXUHV7KHJHRORJLFDOK\GURVWUDWLJUDSKLFPRGHOZDV mainly constructed by A.K. Åberg in co-operation with S.C. Åberg. The manuscript was jointly written by S.C. Åberg and K. Korkka-Niemi, with contributions and comments from the co-authors.

Publications wrote within the degree but not included in thesis

Åberg, A. K., Salonen, V-P., Korkka-Niemi, K., Rautio, A., Koivisto, E. & Åberg, S. C.

2017. GIS-based 3D sedimentary model for visualizing complex glacial deposition in Ker- VLO|)LQQLVK/DSODQG%RUHDO(QYLURQPHQW5HVHDUFKS௅

Åberg, S., Åberg, A., Korkka-Niemi, K., Salonen, V.-P. 2017. Hydrostratigraphy and 3D 0RGHOOLQJRID%DQN6WRUDJH$൵HFWHG$TXLIHULQD0LQHUDO([SORUDWLRQ$UHDLQ6RGDQN\Ol Northern Finland. – In: Wolkersdorfer, C.; Sartz, L.; Sillanpää, M. & Häkkinen, A.: Mine Water & Circular Economy (Vol I). – p. 237–244.

Abbreviations

CLGB Central Lapland Greenstone Belt

DEM Digital elevation model

e.g. exempli gratia

EIA Environmental impact assessment

ELY centre Elinkeino-, liikenne-, ja ympäristökeskus (Centre for Economic Development, Transport and the Environment

EMR Episodic master recession curve method

et al. et alia

GDE Groundwater-dependent ecosystem

GIE Groundwater-iQÀXHQFHGecosystems

GPR Ground penetrating radar

GUI Graphical user interface

K Hydraulic conductivity

LiDAR Light detecting and ranging

m a.s.l. Metres above sea level

NLSF National Land Survey of Finland

SWE Snow water equivalent

UAV-TIR Unmanned aerial vehicle thermal infrared

UPW Upstream Weighting package in MODFLOW

WTF Water table fluctuation method

(14)

/LVWRI¿JXUHV

Fig 1 Location of the study area of this study, mining projects, exploration areas and reservations, page 18

Fig 2 6FKHPDWLFGLDJUDPSUHVHQWLQJWKHZRUNÀRZDQGUHODWLRQVRIWKHPHWKRGV, page 22

Fig 3 3D view of the geological model used in Paper III, page 23

Fig 4 6LPSOL¿HGFRQFHSWXDOPRGHORIJURXQGZDWHUÀRZSDWKVLQULYHUEDQNVDQGLQWKH western part of the mire, page 28

Fig 5 The iterative model development chart of this studySDJHV௅

Fig 6 Cross-sections of the HSM models with hydraulic head contours of four ground ZDWHUÀRZPRGHOVSUHVHQWHGLQ3DSHU,,, page 33

Errors and corrections in published Paper I

In page 160 line 36: m s^-1ĺPV^-2

In page 168 line 19: mm month^-1ĺPP\HDU^-1

(15)

1. Introduction

*URXQGZDWHUÀRZDQGJURXQGZDWHUUHFKDUJH discharge patterns are important parts of the hydrological cycle and need to be acknowledged at mining development sites. Moreover, groundwater discharge areas are essential habitats of groundwater-dependent (GDE) or groundwater- LQÀXHQFHGHFRV\VWHPV*,(DQGWKHLULGHQWL-

¿FDWLRQ LV LPSRUWDQW WR UHGXFH SRWHQWLDO ULVNV and impacts of mining activity on the environ- PHQW0LQLQJDFWLYLWLHVLQÀXHQFHWKHTXDQWLW\

and quality of the water within the mine area and its surroundings, changing the hydrologi- FDODQGWRSRJUDSKLFDOFRQGLWLRQVDQGD൵HFWLQJ JURXQGZDWHUÀRZ3XQNNLQHQet al. 2016). On the other hand, understanding the water balance is critical in planning, designing and positioning mining activities, and this needs to be carefully considered in environmental impact assessment (EIA) (Marandi et al. 2014, Salonen et al. 2014, Krogerus and Pasanen 2016, Punkkinen et al.

2016). In recent years, sustainable development and the management of groundwater resources in PLQLQJD൵HFWHGDUHDVKDYHEHHQDFNQRZOHGJHG (Raghavendra and Deka 2015). Green mining projects should consider environmental aspects in all project phases, from early exploration to mine closure (Nurmi and Wiklund 2012, Nurmi 2017). Attention has also been directed to understanding the geological complexity of mining environments from a hydrogeological perspective (Artimo et al. 2004, Wycisk et al. 2009, Salonen et al. 2014). Understanding the circulation of groundwater and its connections to surface water, wetlands and possible fractured and weathered bedrock aquifers in baseline studies increases the transparency and thus the social licence to oper- ate in environmentally sensitive areas.

The Central Lapland Greenstone Belt

&/*%LVWKHODUJHVWPD¿FYROFDQLFSURYLGHQFH

in Fennoscandia (Eilu et al. 2012) with a high ore potential and considerable mining interest.

Multiple active mines, such as Kevitsa and Kit- tilä, ongoing mining projects such as Sakatti, and numerous claims, exploration permits and reservations are located within the area (Fig. 1). The CLGB coincides with an ice divide area having had relatively weak glacial erosion during last ice ages (Hirvas 1991, Pulkkinen 1983, Ebert et al. 2015). Due to weak glacial erosion, weath- ered bedrock has been preserved on the top of the fractured upper part of the crystalline bedrock (Hirvas 1991, Hall et al. 2015), enabling interaction of the shallow groundwater with ground- ZDWHULQWKHEHGURFN7KH4XDWHUQDU\VXU¿FLDO deposits consist of scattered interlayered till and sorted units, which are typical of river valleys in central Lapland (Lahermo 1970). The combina- WLRQRIVFDWWHUHGÀXYLDODQGJODFLJHQLFVHGLPHQWV with weathered and fractured bedrock creates complex aquifer-aquitard systems in extensive mine prospecting and reservation areas in the

&/*%FDXVLQJJURXQGZDWHUÀRZSDWWHUQVWKDW are challenging to estimate. The complexity of shallow groundwater systems and their connection to weathered and fractured bedrock should be carefully considered in the environmental impact assessment of mining projects located in weak glacial erosion areas.

Surface water bodies are an abundant feature in central Finnish Lapland, and rivers with their outwash plains are common (Lahermo 1973).

Operating mines and exploration areas (Fig. 1) are often located close to surface water bodies and open mire ecosystems. However, the rivers in Lapland are rarely in a natural state. The construction of the hydroelectric power plants and regulation of the rivers have changed the hydrology of the catchment areas. Moreover, most of the mires in Finland have been drained since the 1950s (Aapala 2001, Ruuhijärvi and Lind- KROPZKLFKKDVEHHQREVHUYHGWRD൵HFW

(16)

the groundwater discharge patterns in mire areas (Rossi et al. 2012).

Recognizing the potential habitats of GDE that are commonly related to groundwater–surface water interaction areas or mire areas has become an objective required by the Groundwa- ter Directive (EC 2006) and by law in Finland since 2014 (Act on Water Resources Manage- ment 1263/2014). Groundwater supports a high level of biodiversity and provides habitats in many endangered GDEs in wetlands (e.g., Kløve et al. 2011, Aldous and Bach 2014). The hydrol- ogy of patterned fen/mire complexes depends on VXUIDFHZDWHUUXQR൵VQRZPHOWSUHFLSLWDWLRQ JURXQGZDWHUGLVFKDUJHDQGVSULQJÀRRGV/DS- palainen 1970, Siegel and Glaser 1987, De Mars et al. 1997, Ruuhijärvi and Lindholm 2006, Acre- man and Holden 2013, Isokangas et al. 2017).

6XUIDFHZDWHULQÀRZDQGJURXQGZDWHUGLVFKDUJH enhances the nutrient level in fens, enabling more diverse habitats for plants (Malmer 1986, Sie- gel and Glaser 1987). However, Laitinen et al.

(2005) noted that the changes in groundwater ÀRZSDWWHUQVDQGGLVFKDUJHLQZHWODQGDUHDVVXFK as in fen/mire complexes have been inadequately studied. Furthermore, the composition of the subpeat sediments is important, since the existence of a high conductivity zone beneath a mire en- KDQFHVWKHYHUWLFDOÀRZFUHDWLQJVSDWLDOYDULD- tion in groundwater discharge and recharge areas (Siegel 1988, Reeve et al. 2000). The hydrology RIDPLUHFRPSOH[FDQEHD൵HFWHGE\YDULDWLRQ LQWKHULYHUVWDJHUHODWHGWRQDWXUDOÀRRGLQJRU river regulation if the mire is located in the prox- imity of a river.

*URXQGZDWHUÀRZPRGHOOLQJLVFRPPRQO\

used to study the groundwater interactions with surface water, as well as groundwater recharge, GLVFKDUJHDQGÀRZSDWWHUQV)UHH]HDQG:KLWK- erspoon 1967, Scibek et al. 2007, Barthel and

%DQ]KDI*URXQGZDWHUÀRZPRGHOOLQJLV suitable for sensitive areas, as it can be used to

investigate the water balance and groundwater ÀRZSDWWHUQVDQGPDNHSUHGLFWLRQVIRUIXWXUH changes (Marandi et al. 2013) in a non-destruc- WLYHPDQQHU,QDGGLWLRQJURXQGZDWHUÀRZPRG- el results can reveal the potential environmental risks of anthropogenic activities if the hydrostratigraphic details of the model are adequate.

This hydrostratigraphic detail depends on the research questions, data and time available for model construction (Ross et al. 2005). Modelling codes that can create hydraulic conductivity (K)

¿HOGVEDVHGRQSVHXGRJHQHWLFPRGHOVVXFKDV braided river bar structures can be used to construct detailed models (Bennett et al. 2019). In contrast, simple models are faster to produce and better serve strict schedules. However, if more time is used for the hydrostratigraphic details of a model to increase the knowledge of model area and to reduce the potential future risks and environmental impacts.

The hydrostratigraphic complexity, variation and distribution of hydraulic conductivity zones and topographical variations are the most impor- WDQWIHDWXUHVD൵HFWLQJWKHORFDWLRQRIJURXQGZD- ter discharge and recharge patterns (Freeze and Whitherspoon 1967, Freeze and Cherry 1979, +D\DVKLDQG5RVHQEHUU\,QÀDWO\LQJDU- eas, even small topographic variation can impact on groundwater recharge/discharge patterns (Freeze and Whitherspoon 1967, Reeve et al.

2000). For example, in patterned fens, varia- WLRQLQVWULQJVDQGÀDUNVJHQHUDWHVVPDOOVFDOH groundwater recharge and discharge areas (van der Ploeg et al. 2012). Subsurface complexity af- IHFWVWKHJURXQGZDWHUÀRZSDWWHUQGXHWRUHIUDF- WLRQRIWKHÀRZOLQHVLQORFDWLRQVRIK change (Freeze and Cherry 1979). The refraction causes FRPSOH[ÀRZSDWWHUQVLQKLJKO\YDULDEOHK areas and disperses recharge and discharge areas (Freeze and Cherry 1979, Reeve et al. 2000).

(17)

1.1. Research aims

7KH ¿UVW REMHFWLYH RI WKLV GRFWRUDO VWXG\ ZDV WR$XQGHUVWDQGWKHIDFWRUVD൵HFWLQJUHFKDUJH GLVFKDUJHDQGÀRZSDWWHUQVRIVKDOORZJURXQG- water and bedrock groundwater and their interactions with surface water in the study area in FHQWUDO)LQQLVK/DSODQG%H[DPLQHWKHH൵HFW of hydrostratigraphic complexity on groundwa- WHUUHFKDUJHGLVFKDUJHDQGÀRZSDWWHUQVDQG&

LQYHVWLJDWHWKHH൵HFWRIULYHUUHJXODWLRQRQWKH hydrology of mires and habitats of groundwater- LQÀXHQFHGPLUHVSHFLHVE\FRPSDULQJWKHSUHV- ent situation with the pre-regulation situation. 3D JURXQGZDWHUÀRZPRGHOOLQJDQGÀRRGPRGHO- ling were used to examine the present and past D൵HFWVRQWKHK\GURORJ\RIWKH9LLDQNLDDSDPLUH

The second objective of the study was to FUHDWH D ZRUNÀRZ IRU K\GURVWUDWLJUDSKLF DQG JURXQGZDWHUÀRZPRGHOOLQJDGHTXDWHIRUEDVH- line studies of mining development sites in re- cently glaciated areas with a complex hydrostratigraphy consisting of shallow groundwater interacting with surface water and bedrock groundwater in weathered/fractured bedrock.

2. Description of the study site

2.1. Geological and hydrological background

The study site is located in northern Finland in the municipality of Sodankylä, in the western part of the Natura 2000-protected Viiankiaapa mire in a formerly glaciated area that is characterised by weak glacial erosion (Fig. 1b). A prominent Cu- Ni-PGE mineralization called Sakatti has been found underneath the Viiankiaapa mire (Brown- scombe et al. 2015). The Sakatti deposit is situ- ated in the CLGB, which consists of Paleopro- terozoic sedimentary and volcanic rocks (Brown- scombe et al. 2015), including quartzites, mica schists and gabbros (Tyrväinen 1980, Pulkkinen

1983, Tyrväinen 1983). The uppermost part of the crystalline bedrock is partly highly weathered and fractured (Pulkkinen 1983, Paper II).

According to Lahermo (1973), outwash plains are the most abundant groundwater res- ervoirs of the study area. Weathered and fractured bedrock is covered with alternating Quaterna- ry tills and sorted sediments mainly consisting RIÀXYLDODQGJODFLRÀXYLDOGHSRVLWVDQGDHROLDQ clastic sediments often covered with peat deposits (Åberg A.K. et al. 2017). Groundwater res- HUYRLUVDUHVFDWWHUHGDQGFRQ¿QHGDQGSHUFKHG units also occasionally occur (Paper I). The multiple till layers originate from separate Weichse- lian glacial events (Åberg A.K. et al. 2017). The aquifer-aquitard system can be complex due to existing subtill sands, gravels and gyttja layers, as described in Åberg A.K. et al. (2017). Minero- genic Quaternary deposits are overlain by peat deposits in the Viiankiaapa and Vanttioaapa mire DUHDV7KHDUHDLVWRSRJUDSKLFDOO\UHODWLYHO\ÀDW with altitude variation from about 180 to 207 metres above sea level (m a.s.l.).

The study area is located in a boreal and subarctic region. The mean annual temperature at Sodankylä is -0.4 °C and the June and January average temperatures are +15 °C and -12 °C, UHVSHFWLYHO\ KWWSLOPDWLHWHHQODLWRV¿YXRVLWL- lastot [In Finnish]). The average annual precipitation in the Sodankylä area is 500–650 mm yr^-1, of which about half falls as snow, usually from November to May (http://ilmatieteenlaitos.

¿YXRVLWLODVWRW>,Q)LQQLVK@0HDQHYDSRUDWLRQ was about 320 mm yr^-1 between 1960 and 2011 (Moroizumi et al. 2014).

The water balance of the study area is posi- tive, since annual precipitation exceeds the evaporation rate, and groundwater recharge dom- inates within the valley of the River Kitinen.

The River Kitinen and other smaller rivers are mainly gaining rivers. The main groundwater recharge occurs from April to September. The high-

(18)

190 210

200 230

260 270

180

240

190 190

180

190 210

190

190 190 200

200

27°0'E 26°55'E

26°50'E 26°45'E

26°40'E 26°35'E

67°36'N 67°36'N

67°34'N 67°34'N

67°32'N 67°32'N

67°30'N 67°30'N

67°28'N 67°28'N

Viiankiaapa Kitinen

Sattanen

0 2 4km

Matarakoski hydroelectric power plant

Kersilö

Sakattilampi Viiankijärvi Käppälä- aapa

Särkivaara

Siurunmaa

Sodankylä

Kelukoski hydroelectric power plant

220

200 Kitinen

Vanttioaapa

Moskuvaara

Tiukuoja Sakatti deposit

Model area Paper I Model area Paper II & III Catchment area of Kersilö (65.821) Natura 2000

Contours Water body Gravel intake area Open mire Forested mire

BB

B

Sodankylä

Rovaniemi

32°E 30°E

28°E 26°E

24°E 22°E

20°E

70°N69°N68°N67°N66°N

30° E 30° E

25° E 25° E

20° E 20° E

65° N60° N

Legend

Study area Water body

Catchment of river Kemijoki BActive mine

BUnder development Exploration permission/claim Reservation

Week glacial erosion

± ±

0 100 200km

0 50 100km

Russia Norway

Sweden Norway

Sweden

Russia

Finland

Kemijoki

Lokka reservoir Porttipahta

reservoir

Lake Kemijärvi Kitinen

(c) (b)

(a)

±

(d)

Kevitsa Kittilä

Lampivaara

Kemi Kalkkimaa

Ristimaa Hannukainen

Sokli

Ahmavaara Konttijärvi

Gulf of Bothnia

Figure 1. (a) The location of the study area presented with the catchment area of the Kemijoki river and ongoing PLQLQJSURMHFWVDQGUHVHUYDWLRQV௅E:HDNJODFLDOHURVLRQDUHDVLQGLFDWHDUHDVZLWKSRVVLEO\

complex stratigraphy and preservation of weathered bedrock. (c) The model areas and the study area in closer view.

(d) Outcrop of the Kärväsniemi stratigraphic site (Åberg A.K et al.ZKLFKZDVXVHGDVEDVLVIRUWKHGH¿QLWLRQRI K\GURVWUDWLJUDSKLFXQLWVRIWKHPRGHOV:DWHUERGLHVDQGFDWFKPHQWDUHDVDUHPRGL¿HGDIWHUWKH)LQQLVK(QYLURQPHQW

(19)

,QVWLWXWH7KHDGPLQLVWUDWLYHERUGHUVRI)LQODQGEDVHPDSVDQGWHUUDLQHOHPHQWVZHUHPRGL¿HGDIWHUWKH1DWLRQDO/DQG 6XUYH\RI)LQODQG1/6)0LQHVDQGH[SORUDWLRQDUHPRGL¿HGDIWHUWKH*HRORJLFDO6XUYH\RI)LQODQG([SORUDWLRQ SHUPLVVLRQFODLPVDQGUHVHUYDWLRQV)LQQLVK6DIHW\DQG&KHPLFDO$JHQF\78.(6:HDNJODFLDOHURVLRQLVPRGL¿HG after Ebert et al.6DNDWWLGHSRVLWLVPRGL¿HGDIWHU3|\U\$GPLQLVWUDWLYHERUGHUVRI6ZHGHQ1RUZD\DQG 5XVVLDVRXUFH*$'0)LJXUHGSKRWRJUDSK\$QQLNDcEHUJ7KHFRRUGLQDWHV\VWHPLV(85()),170

est recharge rate in spring results from the spring thaw and the secondary maximum results from intense rainfall in summer and autumn (Paper I).

2.2. Viiankiaapa mire

The Viiankiaapa mire is a nutrient-rich aapa-type SDWWHUQHGIHQZLWKGLVWLQFWLYHVWULQJDQGÀDUNSDW- WHUQVWKDWKDYHIRUPHGSHUSHQGLFXODUWRWKHÀRZ direction of the mire water (Foster et al. 1983).

The Viiankiaapa mire started to develop after the Late Weichelian, about 10 000 years ago (Suon- perä 2016). The sediments underlying the peat mainly consist of aeolian sands and braided river sediments (Lappalainen 2004, Paper I). Sandy sediments are covered with gyttja (Lappalainen 2004) deposited during the Ancylus Lake phase about 10 300–10 200 years ago (Åberg A.K et al.

2017). Gyttja is covered by coarse detritus gyttja deposited before overgrowth of the lake stage (Lappalainen 2004). The lowermost peat layers consist of EU\DOHV peat, which changes upwards to a more mixed composition of carex–EU\DOHV peat with varying proportions of carex (Lappa- lainen 1970). The deposition rate of peat has been about 0.6 mm year^-1 (Lappalainen 2004).

The Viiankiaapa mire forms a habitat for mul- WLSOHSRVVLEO\JURXQGZDWHULQÀXHQFHGSODQWVSH- cies (Hyvärinen et al.,8&1FODVVL¿FDWLRQ in brackets), such as Hamatocaulis vernicosus (near threatened, NT), Hamatocaulis lapponi- cus (vulnerable, VU), Meesia longiseta (endan- gered, EN) and Saxifraga hirculus (vulnerable, VU) (Kulmala 2005, Eurola and Huttunen 2006).

2.3. The River Kitinen

The River Kitinen is a tributary of the Kemijoki river, which discharges into the Gulf of Bothnia (Fig. 1). The River Kitinen has been regulated since the 1970s by hydroelectric power plants

(Alanne et al. 2014), which has changed the hy- GURORJLFDOVHWWLQJVRIWKHIUHHÀRZLQJULYHULQWR regulated pools separated by hydroelectric power plants. The hydroelectric power plants have been used for electricity production and simul- WDQHRXVO\WRUHJXODWHWKHVSULQJÀRRGVWKDWXVHG to cause economic damage to villages along the river shores. Regulation reduced the gradient of the river and simultaneously caused a rise in the river stage behind the dams of the hydroelec- WULFSRZHUSODQWV5HJXODWLRQDOVRD൵HFWHGWKH groundwater–surface water interactions, since the rise in the river stage changed the locations of the groundwater discharge areas in the river EDQNVDQGGHFUHDVHGULYHULQHÀRRGLQJ

3. Materials and methods

3.1. Materials

The study was performed with complex multiple VRXUFHGDWDVHWVWKDWLQFOXGHGGL൵HUHQWREVHUYD- tions from variable sources (Fig. 2), mainly from AA Sakatti Mining Oy, NLSF, the Geological Survey of Finland, Kemijoki Oy, the ELY Centre and several unpublished reports. The quality of the datasets varied considerably, which was considered in model construction. In addition, new data were collected during the study, including groundwater table observations and monitoring GDWDJURXQGSHQHWUDWLQJUDGDU*35SUR¿OLQJ thermal infrared imaging, the stable isotope composition of water and hydrogeochemical data.

The datasets used for constructing the 3D geological models presented in Åberg A.K. et al.

(2017) and Åberg A.K. et al. (2020) were com- piling data from multiple sources, including a till geochemistry dataset, multiple drilling datasets from AA Sakatti Mining Oy, peat drillings, bedrock outcrops, excavated exposures, GPR pro-

(20)

¿OHVDVZHOODVEHGURFNPDSVDQGVXU¿FLDOGH- posit maps from the Geological Survey of Fin- land. The model surface was based on a two- metre-resolution LiDAR DEM (light detection and ranging digital elevation model) produced by NLSF. In addition, the pre-regulation topog- UDSK\ZDVUHFUHDWHGIRUÀRRGPRGHOOLQJDQGIRU JURXQGZDWHUÀRZPRGHOOLQJRIWKHSUHUHJXOD- WLRQJURXQGZDWHUÀRZ$EDVHPDSRI1/6)IURP 1989 and a river cross-section from Huokuna (1991) were used to construct the topography of the gravel intake area and shore area, which was covered with water after regulation (details presented in Paper III).

Groundwater table measurements were available from both automatic monitoring stations and manual observations. The latter included data from reports in 1988 (Lapin vesi- ja ympäristö- piiri 1988, unpublished report) and 1995 (Lapin ympäristökeskus 1998, unpublished report), and PRUHUHFHQWGDWDIURPD¿HOGFDPSDLJQLQ presented in Paper I. Automated monitoring sta- tion data originated from AA Sakatti Mining Oy GDWDVEDVHLQ௅*URXQGZDWHUWDEOHRE- servations were used to calculate and estimate groundwater recharge and also to calibrate all WKH'JURXQGZDWHUÀRZPRGHOV6XUIDFHZDWHU observations were derived from a 1989 map and D/L'$5'(0IURP௅6XUIDFHZDWHU observations from the LiDAR DEM were used WRGH¿QHWKHVXUIDFHZDWHUÀRZGLUHFWLRQVDQG FDOLEUDWHWKHJURXQGZDWHUÀRZPRGHOV

The hydraulic conductivities of Quaternary VHGLPHQWVZHUHGH¿QHGRQWKHEDVLVRIJUDLQVL]H analyses conducted with the Saurbrey method (Sauerbrey 1932) and slug-tests (Golder Associ- ates 2012, unpublished data; SRK, unpublished data, 2019; AA Sakatti Mining Oy, unpublished data, 2019). The hydraulic conductivity of the fractured bedrock was interpolated from hydrogeological testing in boreholes, including packer, spinner and constant head tests (SRK, un-

published data, 2019), along with two slug tests (Golder Associates, unpublished data, 2012).

The locations of modelled groundwater discharge areas in the eastern river shore were con-

¿UPHGZLWKWKHUPDOLQIUDUHGLPDJLQJ7,5DQG the stable isotopic composition of water (Paper ,7,5GDWDZHUHDFTXLUHGGXULQJ¿HOGZRUNZLWK an unmanned aerial vehicle (UAV). Samples for stable isotopic composition analysis were col- OHFWHGLQPXOWLSOH¿HOGZRUNFDPSDLJQVDQGWKH data are presented in Lahtinen (2017), Paper I and in Bigler (2019).

River stage and river discharge data acquired from Kemijoki Oy were used to study the pre- regulation and present hydrological settings of WKHVWXG\DUHDDQGLQÀRRGPRGHOOLQJWKH\ZHUH XVHGWRGH¿QHWKHSUHUHJXODWLRQÀRRGGLVWULEX- tion (Paper III).

The locations of the examined groundwa- WHULQÀXHQFHGVSHFLHVHamatocaulis lapponi- cus, Hamatocaulis vernicosus and Saxifraga hirculus DQG ÀRRGGHSHQGHQW VSHFLHV Carex microglochin, 0RHKULQJLDODWHULÀRUD and Ely- PXVPXWDELOLV) were acquired from the Eliölajit database of the Finnish Environment Institute (SYKE, unpublished data) and the AASMOy biological database (2021, unpublished data).

The distributions of the aforementioned species were compared with the modelled regula- WLRQLQÀXHQFHGFKDQJHLQWKHJURXQGZDWHUWDEOH changes in groundwater discharge areas between pre- and post-regulation models, and the simu- ODWHGSUHUHJXODWLRQÀRRGFRYHUDJH3DSHU,,, Simultaneously, the relationship between the altitude of the groundwater table and the present distribution of the plant species was examined.

3.2. Methods

This doctoral study was completed with multiple GL൵HUHQWPHWKRGVWKDWZHUHXVHGWRFKDUDFWHU- ize the present and pre-regulation hydrological settings of the study area and evaluate its vul-

(21)

nerability to mining practices. The main methods included 1) hydrostratigraphic 3D modelling 3DSHUV,௅,,,'JURXQGZDWHUÀRZPRGHO- OLQJ3DSHUV,௅,,,XVHGWRGH¿QHDQDGHTXDWH level of hydrostratigraphic detail (Paper II) and WKHSUHVHQWDQGSUHUHJXODWLRQJURXQGZDWHUÀRZ and recharge/discharge patterns (Paper III), 3), ÀRRGPRGHOOLQJ3DSHU,,,DQGFKDUDFWHUL]D- WLRQRIWKHK\GURORJLFDOIHDWXUHVWKDWPD\D൵HFW the threatened GIE within the study area (Paper ,,,0DMRUZRUNÀRZVDUHSUHVHQWHGLQ)LJ 3.3. 3D Hydrostratigraphic

modelling with Leapfrog Geo (Papers I, II and III)

The hydrostratigraphic models based on 3D geological modelling were generated with Leapfrog Geo (Seequent Ltd 2020). The modelling, based on a combination of explicit modelling and im- plicit modelling, is presented in more detail in Åberg et al. (2017) and Paper II. The datasets used to generate the hydrostratigraphic units were mainly from the Kersilö database (Åberg A.K. et al. 2017) and drilling data and hydraulic conductivity data were obtained from AA Sakatti Mining Oy. The 3D geological models were sim- SOL¿HGDVK\GURVWUDWLJUDSKLFPRGHOVE\H[FOXG- ing units with the smallest and most uncertain GLVWULEXWLRQV3DSHUV,௅,,,7KH'JHRORJLFDO model was updated multiple times based on the UHVXOWV RI JURXQGZDWHU ÀRZ PRGHOV QHZ RE- servations and drilling data that were acquired during this study.

The selection of hydrostratigraphic units for WKHVXU¿FLDOGHSRVLWVZDVPDLQO\EDVHGRQWKH separation of glacial till (3 units) from interlayered sorted sediments (4 units). In addition, a basal unknown sediment unit and two peat units were used as hydrostratigraphic units (Fig. 3).

It was assumed that the hydraulic conductivity of till is generally lower than that of sorted sediments (1.2 × 10^-6–4.3 × 10^-3 m s^-1; Paper

II). However, the K of till was highly variable (5.0×10^-8 to 8.3×10^-5 m s^-1; Paper II), and some of observed high conductivity sandy tills act as aquifers rather than aquitards. Simultaneously, it was acknowledged that major hydrostratigraphic units consisted of heterogeneous mate- rial, which was considered in the parametrization RIWKHJURXQGZDWHUÀRZPRGHOV7KHEHGURFN ZDV¿UVWFUHDWHGZLWKRQH3DSHU,WKHQZLWK two (Paper II) and lastly with three units (Paper III). The bedrock in Paper III is comprised of clay-type weathered bedrock with assumed low conductivity, grus-type weathered bedrock with assumed intermediate hydraulic conductivity and a fractured bedrock unit with variable hydraulic conductivity based on packer and spinner tests.

The fractured bedrock unit also included separate parameter units for major discrete faults FODVVL¿HGLQWRWKUHHJURXSVWKUXVWIDXOWVEULWWOH IDXOWVDQGFOD\¿OOHGIDXOWV7KHDVVXPSWLRQZDV that brittle faults have a relatively high K, thrust faults have a lower KDQGFOD\¿OOHGIDXOWVWKH lowest K, which corresponded to the K of a clay weathered bedrock unit (Paper III). The peat of the mire was assumed to be homogeneous for VLPSOLFLW\LQWKH¿UVWK\GURVWUDWLJUDSKLFPRGHO presented in Paper I. In later models (Papers II and III), the mire was modelled with two layers, including an upper less decomposed and lower more decomposed layer, edited in the groundwa- WHUÀRZPRGHOJULGLQWKH0RGHO0XVHJUDSKLFDO user interface (GUI) (Winston 2009). The hydrostratigraphic models were revised multiple times EDVHGRQQHZREVHUYDWLRQV*35SUR¿OLQJGULOO- LQJGDWDDQGJURXQGZDWHUÀRZPRGHOUHVXOWV 'JURXQGZDWHUÀRZ

modelling with MODFLOW- NWT (Papers I, II and III)

The MODFLOW series (McDonald and Har- EDXJK LV D ¿QLWH GL൵HUHQFH PRGHOOLQJ FRGHFDOFXODWLQJWKHJURXQGZDWHUWDEOHDQGÀRZ

(22)

3DSHU, 3DSHU,, 3DSHU,,,

*URXQGZDWHUIORZPRGHOOLQJ 02')/2:1:7

Flood modelling HECRAS Input data in

3DSHU,

• *HRORJLFDO model 3.5 x 3 km

• Hydraulic conductivities

• *URXQGZDWHU recharge

• Boundary conditions

• River stage

• Fractured bedrock

,QSXWGDWD

• 3UHregulation topography

• River discharge data from 1966

• River stage data 1987

• 3UHFLSLWDWLRQ DQG6:(IURP 1966 ad 1987 Additional input data

LQ3DSHU,,

• *HRORJLFDO models

• *0DQG*0

• Hydrostratigraphic models HSM14

• :HDWKHUHG bedrock unit including gruss type and clay type weathering

Model calibration 3(67

Additional input data in 3DSHU,,,

• 35(DQG

3267PRGHOV

• 3UHregulation topography

• 3UHregulation river stage

• Three layer bedrock (clay, gruss, fractured bedrock)

• *URXQGZDWHUhead observations from 1988

0RGHORXWSXW 3DSHU,

• *URXQGZDWHUWDEOHIORZSDWWHUQV

• Recharge and discharge areas

$GGLWLRQDOLQ3DSHU,,,

• Changes in

groundwater table and recharge/discharge patterns due to regulation

$GGLWLRQDOLQ3DSHU,,

• Effect of model complexity to recharge and discharge patterns

• Hydraulic connections between surficial deposits and weathered/fractured bedrock

0RGHORXWSXW

• flood coverage in 1966 and in 1987

*URXQGZDWHU head observation Crosssections Drillhole data

*URXQGZDWHUKHDGobservations

*URXQGpenetrating radar

Hydraulic conductivity measurements LiDAR DEM

Outcrops, field observations Surficial deposit and bedrock maps

Relations to groundwater influenced and flooddependent ecosystems Model

verification with stable isotopes of water and UAVTIR images 6 x 8 km models

Sensitivity analysis

Figure 2.6FKHPDWLFGLDJUDPSUHVHQWLQJWKHZRUNÀRZDQGUHODWLRQVRIWKHPHWKRGVXVHGLQ3DSHUV,,,DQG,,,

in each discrete cell, which has been used for

\HDUVLQJURXQGZDWHUÀRZPRGHOOLQJ7KH02'-

FLOW-NWT version (Niswonger et al. 2011), which uses the Upstream Weighting (UPW)

(23)

package with the Newton formulation, was cho- sen as the modelling code, since it enables the PRGHOOLQJRIXQFRQ¿QHGDTXLIHUVZLWKDKLJK- ly variable layer thickness and scattered uncon-

¿QHGXQLWV7KH83:SDFNDJHSUHYHQWVDFHOO from being fully dry and keeps all cells active, which enables water table variation in a vertical direction between two or more stacked cells during model iterations. In contrast, in the widely used MODFLOW-2005 (Harbaugh 2005), wet- ting and drying of the cells is a challenging is- sue in thin scattered units. This is problematic, HVSHFLDOO\LQWKLQOD\HUHGXQFRQ¿QHGV\VWHPV which are common in formerly glaciated areas.

MODFLOW-NWT has been successfully used with thin aquifers (Hunt and Feinstein 2012), and the results have been compared with analytical solutions (Zaidel 2013).

The hydrostratigraphic models in Papers I–

,,, ZHUH ¿UVW FRQVWUXFWHG ZLWK /HDSIURJ *HR (Seequent 2020 Ltd) and then converted to MODFLOW grids with the Hydrogeology tool of Leapfrog (further details in Papers I and II).

The parameter zones of major hydrostratigraphic units were created in Leapfrog, and additional parameter zones were added in the ModelMuse GUI (Winston 2009), which was used for further JURXQGZDWHUÀRZPRGHOFRQVWUXFWLRQ3DSHUV,௅

5LYHU.LWLQHQ 3HDW 7RSGHSRVLWV 8SSHUWLOO

0LGGOHVRUWHGGHSRVLWV 0LGGOHWLOO

/RZHUWLOODVVHPEODJH /RZHVWVDQGVDQGJUDYHOV

%DVDOXQNQRZQVHGLPHQW

&OD\W\SHZHDWKHUHGEHGURFN

*UXVW\SHZHDWKHUHGEHGURFN )UDFWXUHGEHGURFN

/HJHQG

/RZHUVRUWHGGHSRVLWV

D E

PDVO PDVO

NP NP

9HUWLFDOH[DJJHUDWLRQ[

1 1

Figure 3. *HRORJLFDOPRGHORIWKHVWXG\DUHDIURP3DSHU,,,FRQVWUXFWHGZLWK/HDSIURJ*HRLQDH[SORGHGYLHZ and (b) 3D view.

III). The hydraulic conductivities of the param- HWHU]RQHVZHUHGH¿QHGZLWKLQWHUSRODWLRQDQG estimation from the results of grain-size analyses conducted with the Sauerbrey method (Sau- erbrey 1932), slug tests, packer tests and spinner tests (further details in Papers I–III). Groundwa- ter recharge was estimated based on the results RIZDWHUWDEOHÀXFWXDWLRQ:7)0HLQ]HU and episodic master recession curve (EMR) (Nimmo et al. 2015) methods.

7KHJURXQGZDWHUÀRZPRGHOVLQ3DSHUV,±

III were calibrated with PEST (Doherty 2015).

The parameter estimation method was used in Papers I and II and the stochastic random parameter method RANDPAR1 (Doherty 2018) in Paper III. In parameter estimation, PEST was XVHGWR¿QGWKHVLQJOHRSWLPDOSDUDPHWHUVHWWKDW KDGWKHEHVWSRVVLEOH¿WEHWZHHQREVHUYHGDQG simulated groundwater heads. In contrast, in the random parameter method RANDPAR1, multi- SOHPRGHOVZLWKGL൵HUHQWSDUDPHWHUVDUHUXQDQG WKHLUUHVXOWVDUHLQVSHFWHGWR¿QGWKHEHVWPRGHO ZKLFKKDVDVJRRGD¿WEHWZHHQREVHUYDWLRQV and simulated heads as possible and reasonable parameter values compared to the measured or estimated K.

7KHJURXQGZDWHUÀRZPRGHOUHVXOWVZHUH YHUL¿HGZLWKWKHVWDEOHLVRWRSLFFRPSRVLWLRQRI

(24)

water (Papers I and III) and with particle tracking (Paper II). In Paper II, MODPATH (Pollock ZDVXVHGWRYHULI\WKHÀRZSDWKVRIWKHRE- served isotopic compositions of water. The sim- XODWHGÀRZSDWKVZHUHFRPSDUHGZLWKWKHLVR- topic composition of springs on the river shore.

The simulated groundwater discharge areas were also visually compared with the distribution of JURXQGZDWHULQÀXHQFHGVSHFLHVLQWKHPLUHDU- ea, and with groundwater discharge areas on the river shore observed with UAV-TIR in Paper I.

'ÀRRGPRGHOOLQJZLWK +(&5$6DQGÀRRGIUHTXHQF\

analysis (Paper III)

Flood modelling was used to study the pre-reg- XODWLRQÀRRGFRYHUDJHRIWKH5LYHU.LWLQHQDQG WRH[DPLQHLWVH൵HFWRQWKHK\GURORJ\RIWKH9LL ankiaapa mire. Flood modelling was performed using HEC-RAS with 2D models (Brunner 1995, Brunner 2016). The HEC-RAS model calculated the river stage based on river discharge data with WKHGL൵XVLRQZDYHHTXDWLRQ%UXQQHU7KH Manning’s roughness value of the river bottom is the most important parameter, which can be used to calibrate the model if river stage or discharge observations are available.

Flood modelling was performed for two ÀRRGV±^th May 1966 and 12–25^th May 1987, based on available river discharge (Kemijoki Oy, unpublished data) and river stage data (Huokuna 1991), respectively. River discharge data from 1966 were used to create a rating curve, which was then used to convert the 1987 river stage data to 1987 river discharge data. The converted 1987 river discharge data were used to simu- ODWHWKHÀRRGFRYHUDJHRI0D\ZKLFKZDV WKHQFRPSDUHGZLWKWKHQDWXUDOVWDWHRIÀRRG- LQJGH¿QHGZLWKÀRRGIUHTXHQF\DQDO\VLV7KH UHFXUUHQFHLQWHUYDORIWKHÀRRGZDVFDO- culated with Log Pearson type III and Weibull methods (Viessman and Lewis 2003). In addi-

tion, the recurrence interval was calculated for the maximum river discharge value of the 1987 ÀRRGWRH[DPLQHLWVUHODWLRQVKLSZLWKSUHUHJX- ODWLRQÀRRGLQJ

3.6. TIR imaging with UAV (Paper I) The TIR survey of groundwater discharge ar- HDVZDVEDVHGRQWKHWHPSHUDWXUHGL൵HUHQFHEH- tween discharging groundwater and surface water (Rautio et al. 2018). Since the groundwater temperature is rather stable and usually close to the annual average air temperature, groundwater discharge can be observed as an anomaly compared to the surface water if the temperature dif- ference is large enough.

The TIR survey was used to identify groundwater discharge areas in the shore of the River Kitinen. The survey was performed with a UAV equipped with TIR and RGB cameras. The tech- nical details are presented in Paper I. The UAV- TIR study was performed in August, when the GL൵HUHQFHEHWZHHQWKHVXUIDFHZDWHU!&

3DSHU , DQG JURXQGZDWHU WHPSHUDWXUHV ௅

°C; AA Sakatti Mining Oy, unpublished data, 2015) was greatest, and groundwater discharge areas were possible to observe as cold anoma- lies from TIR images. The observed groundwater discharge areas were compared with groundwa- WHUÀRZPRGHOUHVXOWVWRHYDOXDWHWKHPRGHOV

3.7. Stable isotopic composition of water (Papers I, II and III)

7KHVWDEOHLVRWRSLFFRPSRVLWLRQVį¹⁸2DQGį' as well as d-excess values, have been used for years to determine the surface water and groundwater fractions from water samples (Hunt et al.

1998, Clay et al. 2004, Rautio 2015). When water is exposed to evaporation, its starts to be depleted of lighter isotopes and its composition deviates from the isotopic composition of meteoric water, which usually has a similar composition to JURXQGZDWHU*DWDQG*RQ¿DQWLQL6XUIDFH

(25)

water bodies are usually exposed to evaporation and fractionated with heavier isotopes. If the iso- WRSHFRPSRVLWLRQHQGPHPEHUVVX൶FLHQWO\GLI- fer, the isotopic composition can be used to dis- tinguish groundwater from surface water (Gat DQG*RQ¿DQWLQL.UDEEHQKRIWet al. 1990, Kendall and Coplen 2001, Rautio and Korkka- Niemi 2015).

7KHVWDEOHLVRWRSHVį¹⁸2DQGį'DQGd-excess YDOXHVZHUHXVHGWRFKDUDFWHUL]HGL൵HUHQWZDWHU types and their possible origin within the study area in Paper I. Stable isotopes were also used WRYHULI\WKHJURXQGZDWHUÀRZPRGHOLQ3DSHUV II and III. The stable isotopic composition was assumed to remain unchanged during groundwa- WHUWUDQVSRUW*DWDQGWKXVUHÀHFWWKHFRP- position of recharging water. Majority of isotope samples used in this study were taken during late summer and early autumn in order to ensure that the isotope signal from snow melt was negligible (Isokangas et al. 2017). In Paper II, MODPATH (Pollock 1994) particle tracking was used with VWDEOHLVRWRSHVWRLQYHVWLJDWHWKHÀRZSDWKVRI GL൵HUHQWLVRWRSHFRPSRVLWLRQVDQGWRFRPSDUH the results with observed stable isotope composition values in groundwater discharge areas in the river banks. In Paper III, stable isotopes of ZDWHUZHUHXVHGWRYHULI\WKHJURXQGZDWHUÀRZ model results in the Matarakoski area (Fig. 1) to examine how the construction of the Matarakoski SRZHUSODQWKDVD൵HFWHGWKHSUHVHQWJURXQGZD- WHUÀRZSDWKV

4. Results

4.1. Results of Paper I

7KHJURXQGZDWHUÀRZPRGHOOLQJUHVXOWVLQGL- FDWHGWKDWPLUHZDWHUIURP9LLDQNLDDSDÀRZV towards the River Kitinen and discharges locally within the mire area and in the shores of the riv- HU7KHJURXQGZDWHUÀRZGLUHFWLRQLQWKHZHVW- ern part of the Viiankiaapa mire was observed

to gradually change during spring from an E–W direction to a more NE–SW direction. Ground- ZDWHUWDEOHÀXFWXDWLRQZDVUHODWHGWRWKHVSULQJ thaw and higher recharge rate in the northern part of the model area. The model results indicated that groundwater recharge in the mire area was related to topographical variation and subpeat sediments. Groundwater recharge was observed in groundwater monitoring wells as groundwater table variation in the river banks and in the mire area.

Complex hydrostratigraphy considerably af- IHFWHGWKHJURXQGZDWHUÀRZSDWWHUQVVLQFHKLJK variation in hydraulic conductivities was present, and multiple scattered till units were observed with scattered sorted sediment units. In addition, bedrock fractures and faults locally af- IHFWHGWKHJURXQGZDWHUWDEOHDQGÀRZDQGGLV- charge patterns.

The stable isotopic composition of water and UAV-TIR survey indicated that groundwater discharge occurred in the eastern shore of the Riv- er Kitinen and eastern fringe of the Viiankiaa- SDPLUH6WDEOHLVRWRSHVDQGJURXQGZDWHUÀRZ modelling results indicated that the heterogene- ity of the fractured bedrock and peat should be studied further in future modelling to improve understanding of the hydraulic connection between the mire and fractured bedrock within the study area.

4.2. Results of Paper II

7KHH൵HFWRIK\GURVWUDWLJUDSKLFFRPSOH[LW\ZDV HYDOXDWHGZLWKIRXUGL൵HUHQWK\GURVWUDWLJUDSKLF models, which indicated that the addition of geo- ORJLFDOGHWDLOWRWKHJURXQGZDWHUÀRZPRGHOLQ- FUHDVHGWKHPRGHO¿WHVSHFLDOO\LIJURXQGZDWHU PRQLWRULQJZHOOVDUHVFUHHQHGDWGL൵HUHQWOLWKR- ORJLFDOXQLWVLQGL൵HUHQWGHSWKV([SOLFLWPRGHO- OLQJZLWK/HDSIURJ*HRZDVIRXQGWREHDQH൵HF- tive but time-consuming method for modelling highly heterogeneous Quaternary systems with

(26)

weathered and fractured bedrock. MODFLOW- NWT was found to be suitable for modelling JURXQGZDWHUÀRZZLWKFRPSOH[4XDWHUQDU\VHG- iment package, since it allowed modelling with XQFRQ¿QHGVFDWWHUHGXQLWVKDYLQJDKLJKO\YDUL- able layer thickness.

The results of the four models indicated that the addition of hydrostratigraphic detail created more complex groundwater recharge and dis- FKDUJHSDWWHUQVDQGLQFUHDVHGWKHYHUWLFDOÀRZ if high variation in hydraulic conductivity was SUHVHQW,QFRQWUDVWJURXQGZDWHUÀRZZDVSUH- dominantly horizontal in models in which Qua- WHUQDU\VHGLPHQWVDQGEHGURFNZHUHVLPSOL¿HGDV one layer per unit. The calibration of the complex models was more challenging, and low sensitive parameters needed to remain unchanged during calibration, causing greater uncertainty in the calibrated parameters and the modelling results.

The creation of complex models was rec- ommended if the research interest is related to GHWDLOHGXQGHUVWDQGLQJRIJURXQGZDWHUÀRZSDW- terns in highly heterogeneous systems, or their interaction and connections with surface waters, and enough time and data are available for the modelling.

4.3. Results of Paper III

The regulation of the River Kitinen and especially the construction of the Matarakoski and Kelukoski hydro-electric power plants in 1995 DQGUHVSHFWLYHO\D൵HFWHGWKHK\GURORJLFDO settings of the Viiankiaapa mire by reducing the hydraulic gradient towards the river. The results RIÀRRGPRGHOVLQGLFDWHGWKDWUHJXODUÀRRGLQJ D൵HFWHGWKHZHVWHUQPRVWSDUWRIWKHPLUHEHIRUH the regulation of the river. The recurrence interval RIWKHÀRRGZLWKDPD[LPXPGLVFKDUJH rate exceeding 900 m³ s^-1, was 11 years, and the UHFXUUHQFHLQWHUYDORIDÀRRGFRPSDUDEOHWRWKDW in 1987, with maximum discharge rate exceeding 600 m³ s^-1, was two years before river regu-

ODWLRQ7KHÀRRGFRYHUDJHPRGHOIRUFRU- UHVSRQGHGZLWKREVHUYHGÀRRGSODLQVHGLPHQWV The rise in the river stage reduced the hydraulic gradient towards the river, increasing the groundwater table in the river banks and western part of the Viiankiaapa mire. River regulation increased the groundwater discharge areas in mire, changed the spring locations in the river shore to a higher altitude, and reduced groundwater discharge into the River Kitinen.

7KHJURXQGZDWHUÀRZPRGHOVDQGWKHVWD- ble isotopic composition of water indicated that the constructed Matarakoski dam changed the JURXQGZDWHUÀRZGLUHFWLRQVRQO\ORFDOO\QH[WWR GDP7KHVWXGLHGÀRRGGHSHQGHGSODQWVSHFLHV GHFOLQHGGXHWRWKHUHGXFWLRQRIÀRRGLQJPD[L- mums after river regulation. The regulation-in- ÀXHQFHGJURXQGZDWHUWDEOHULVHFRYHUHGDOPRVW half of the current habitats of all studied ground- ZDWHULQÀXHQFHGPLUHSODQWVSHFLHV7KHRFFXU- rences of Hamatocaulis vernicosus and Hamato- caulis lapponicus appeared to be related to the high water table and groundwater discharge areas in Viiankiaapa. The model results indicated that WKH9LDDQNLDDSDPLUHLVSDUWO\KXPDQLQÀXHQFHG due to regulation of the River Kitinen.

5. Discussion

5.1. General overview of the JURXQGZDWHUÀRZSDWWHUQV in the study area

Groundwater circulation in the study area is highly complex due to the heterogeneous hydrostratigraphy. Groundwater recharge and discharge within the mire area are controlled by the composition of subpeat sediments and peat, as well as topographical variation. Groundwater re- charges in the river banks and in the forested areas and topographical highs within the mire area.

*URXQGZDWHUPDLQO\ÀRZVLQKLJKFRQGXFWLYLW\

]RQHVRIÀXYLDOVDQGVDQGJUDYHOV)LJDQG

3D hydrogeological modelling at a mining development site with complex hydrostratigraphy in northern Finland