Dynamic model of a municipal wastewater stabilization pond in the arctic

Accepted Manuscript Dynamic model of a municipal wastewater stabilization pond in the arctic Didac Recio-Garrido, Yehuda Kleiner, Andrew Colombo, Boris Tartakovsky PII:

S0043-1354(18)30595-5

DOI:

10.1016/j.watres.2018.07.052

Reference:

WR 13951

To appear in:

Water Research

Received Date: 4 January 2018 Revised Date:

19 July 2018

Accepted Date: 20 July 2018

Please cite this article as: Recio-Garrido, D., Kleiner, Y., Colombo, A., Tartakovsky, B., Dynamic model of a municipal wastewater stabilization pond in the arctic, Water Research (2018), doi: 10.1016/ j.watres.2018.07.052. This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

AC C

EP

TE D

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

ACCEPTED MANUSCRIPT

1

Dynamic Model of a Municipal Wastewater Stabilization Pond in the Arctic

2

4 5 6

Didac Recio-Garrido1, Yehuda Kleiner2, Andrew Colombo2, and Boris Tartakovsky1* National Research Council of Canada 1 6100 Royalmount Ave, Montreal, QC, Canada H4P 2R2 2 1200 Montreal Rd, Ottawa, ON, Canada K1A 0R6

RI PT

3

7

Waste stabilisation ponds (WSPs) are the method of choice for sewage treatment in most arctic

9

communities because they can operate in extreme climate conditions, require a relatively modest

SC

8

investment, are passive and therefore easy and inexpensive to operate and maintain. However,

11

most arctic WSPs are currently limited in their ability to remove carbonaceous biochemical

12

oxygen demand (CBOD), total suspended solids (TSS) and ammonia-nitrogen. An arctic WSP

13

differs from a ‘southern’ WSP in the way it is operated and in the conditions under which it

14

operates. Consequently, existing WSP models cannot be used to gain better understanding of the

15

arctic lagoon performance. This work describes an Arctic-specific WSP model. It accounts for

16

both aerobic and anaerobic degradation pathways of organic materials and considers the periodic

17

nature of WSP operation as well as the partial or complete freeze of the water in the WSP during

18

winter. A uniform, multi-layer (ice, aerobic, anaerobic and sludge) approach was taken in the

19

model development, which simplified and expedited numerical solution of the model, enabling

20

efficient model calibration to available field data.

22 23

TE D

EP

AC C

21

M AN U

10

Keywords: Arctic, facultative lagoon, WSP, dynamic model; icecap

* corresponding author phone: 514-496-2664 e-mail: [email protected]

ACCEPTED MANUSCRIPT 2 24 25

1. INTRODUCTION Waste stabilisation ponds (WSPs, also called “sewage lagoons” or “facultative lagoons”) are used for secondary treatment of municipal sewage by many small communities around the

27

globe because they require a relatively modest investment, they are easy and inexpensive to

28

operate and maintain by locally available personnel. WSPs are the method of choice for sewage

29

treatment in most Canadian arctic communities because they can operate in extreme climate

30

conditions and avoid some of the key challenges of mechanical treatment in the north (Johnson

31

et al. 2014). Despite their advantages, most arctic WSPs are currently limited in their ability to

32

remove carbonaceous biochemical oxygen demand (CBOD), total suspended solids (TSS) and

33

ammonia-nitrogen. WSPs often appear as an immediate upstream process to wetlands, so that

34

direct discharge of sewage to wetlands could occur due to poor WSP performance. Importantly,

35

the wetlands provide 47 - 94% cBOD and 39 – 98% TSS removal (Yates et al. 2012).

SC

M AN U

There are some fundamental differences between WSPs operated in the Arctic and those

TE D

36

RI PT

26

operated elsewhere. Whereas in non-Arctic settings the lagoon is generally full (constant water

38

volume) and is operated as a flow-through system, the extreme sub-zero temperatures in the

39

Arctic winter do not usually allow for any outflow of effluent, as this would freeze instantly.

40

Moreover, as the sewage in the arctic lagoon is ice-capped during most of the year, a low organic

41

removal rate can be expected during this cold period. As a result, the typical Arctic WSP has a

42

volume designed to accommodate an entire year of sewage production with the sewage flowing

43

into the WSP during the year. The WSP is typically emptied (decanted) in late summer to early

44

fall (in some locations WSPs are decanted slowly during the entire summer). After decanting, the

45

WSP is ready to receive the next year’s sewage production. These conditions make all existing

46

WSP models unsuitable for simulating Arctic WSPs.

AC C

EP

37

ACCEPTED MANUSCRIPT 3 Despite the simplicity of WSP design and the practical experience gained from its use for

48

several decades all around the world, a complete understanding of all the physical, chemical and

49

biological processes involved is still an ongoing area of research (Sah et al. 2012). The

50

complexity of the various processes that take place simultaneously in WSP necessitated

51

simplification in modeling efforts. Some existing models focus on specific processes governing

52

lagoon physics/chemistry: for example, many of the existing models focus on hydrodynamics for

53

a better WSP geometry design in 1D, 2D or 3D (Aldana et al. 2005, Martínez et al. 2014, Sah et

54

al. 2011, Salter et al. 2000, Sweeney et al. 2003). Other models take a particular modelling

55

approach, such as completely mixed or plug flow, and focus on the biochemical processes for a

56

better understanding of the biodegradation processes and effluent quality (Beran and Kargi 2005,

57

Dochain et al. 2003, Houweling et al. 2005, Peng et al. 2007). Other models still, focus primarily

58

on describing the sedimentation mechanisms of the particulates (Jupsin and Vasel 2007, Toprak

59

1994) or the daily (Kayombo et al. 2000) or seasonal (Banks et al. 2003, Chaturvedi et al. 2014)

60

oxygen dynamics.

SC

M AN U

TE D

61

RI PT

47

Temperature profiles have also been the object of study and modeling in WSP operated in moderate climates (Gu and Stefan 1995, Sweeney et al. 2005). More recently, the coupling of

63

hydrodynamic equations with biochemical processes produced a complex 3D model (Sah et al.

64

2011) that provides a detailed component distribution in the lagoon. However, such 3D models

65

are difficult to calibrate, given existing scarce experimental data, often lacking measurements at

66

different lagoon depths and locations.

AC C

EP

62

67

Some kinetic and stoichiometric parameters for the biochemical processes can be

68

assumed from well-established models used in other domains, e.g. ADM1 (Batstone et al. 2002)

69

and ASM3 (Henze et al. 2000). However, in general the literature reflects scant information with

ACCEPTED MANUSCRIPT 4 regard to calibration and validation of models with full-scale WSP data (Sah et al. 2012). In the

71

Arctic, this dearth of WSP field data is even greater due to extreme winter conditions and

72

remoteness of WSP locations (Ragush et al. 2015). The formation of a thick ice cover in Arctic

73

lagoons for 6 to 8 months of the year (depending on location) increases modeling complexity as

74

it affects all physical and biological mechanisms of COD removal. The proposed dynamic model

75

endeavours to close the gap by considering the periodic nature of lagoon operation in the Arctic

76

as well as the impact of ice cover on WSP performance.

2. MODEL FORMULATION

79

2.1 Modeling approach

80

SC

78

M AN U

77

RI PT

70

The proposed model uses a multi-layer approach to describe the existence of several distinct zones in a facultative lagoon. The multi-layer approach to modeling was successfully

82

applied in the past to various biological systems with significant oxygen, carbon source, or

83

biomass gradients (Rauch et al. 1999, Tartakovsky and Guiot 1997). Accordingly, the model

84

proposes the existence of three distinct layers at any given time, each assumed to be completely

85

mixed. In the relatively warm period in the summer, the top layer is assumed to be aerobic liquid,

86

the middle layer anaerobic liquid, and the bottom layer is sludge, assumed to always be under

87

anaerobic conditions. The depth of the aerobic layer is expected to be restricted to 20 – 50 cm

88

due to limited oxygen penetration (no forced aeration) and fast oxygen consumption by aerobic

89

bacteria (Hunik et al. 1994, Stephenson et al. 1999, Tartakovsky et al. 2006). The assumption of

90

a thin aerobic layer is also supported by Ragush et al (2015), who observed a consistently low

91

dissolved oxygen concentration throughout the lagoon depth.

AC C

EP

TE D

81

ACCEPTED MANUSCRIPT 5 92

In the cold period the WSP is assumed to be completely ice-capped, which means that there is no oxygen flow from the air into the liquid (the depletion of oxygen in the liquid due to

94

icing is assumed to be nearly instantaneous). Consequently, there is no aerobic layer during the

95

cold season. The anaerobic layer of liquid can become partially or completely frozen, depending

96

on ambient temperatures. The sludge layer is assumed to remain in a liquid state due to its depth

97

as well as due to exothermic microbial activities. Microbial activities such as anaerobic

98

fermentation of organic materials and anaerobic conversion of fermentation products to

99

biomethane are exothermic and result in heat production (Conrad 1999), which helps to prevent

100

the sludge layer from freezing. It is noted however that in very shallow lagoons the sludge could

101

freeze at very low temperatures. In the proposed model such an event can be simulated by setting

102

the growth and biodegradation rates to zero.

SC

M AN U

103

RI PT

93

Figure 1 shows a schematic diagram of the layers considered during the warm and cold periods. During the warm season the depth (or thickness) of the aerobic layer (zAE) is considered

105

to be constant, while the depth of the anaerobic (zAN) and the sludge (zSL) layers increases as a

106

result of the daily addition of raw wastewater.

TE D

104

It is to be noted that settling of solids is not explicitly considered in the proposed model,

108

with the following rationale. The time scale for settling is typically in the range of several hours

109

to a few days at most. Initial attempts at utilizing a much more complex 3D model, which also

110

accounted for settling of solids, showed that at a large time scale (such as one year of lagoon

111

operation), instant settling of the solids can be assumed without significantly affecting modeling

112

results. Also, COD transfer between ice and liquid during ice formation is neglected.

AC C

EP

107

ACCEPTED MANUSCRIPT 6 113 114

2.2 Icecap modeling A separate model was developed to estimate the thickness of ice as a function of the ambient temperature profiles (degree days). This model was developed based on work by Ashton

116

et al (Ashton 1986, 1989) with modifications to account for the heat flux from the relatively

117

warm sewage inflow to the ice layer and for the insulating properties of snow cover. A detailed

118

description of the icecap model is provided in Supplementary Information (Appendix A). The

119

direct coupling of the icecap model to the multi-layer lagoon model would have resulted in more

120

complex computations and a slower combined model with a marginal improvement in accuracy.

121

Moreover, year-round lagoon water temperature data are currently not available but the onset of

122

ice formation and the end of ice thaw periods are easily observed and therefore more often

123

known. The icecap model is therefore used to estimate the time at which the icecap begins to

124

form, the time at which it reaches maximum thickness, and the time at which it completely melts.

125

The maximum icecap thickness is also estimated. These estimated values are subsequently used

126

in the lagoon model under the simplifying assumption of constant rate of freeze (ice growth) and

127

constant rate of ice thaw. It is noted that while the model assumes uniform thickness of the ice

128

layer, at the discharge location of warm sewage to the lagoon the ice layer is either much thinner

129

or nonexistent. However, this localized thinning is very small compared to the entire surface area

130

of the lagoon and is therefore deemed to be negligible.

131

2.3 Biodegradation modeling

SC

M AN U

TE D

EP

AC C

132

RI PT

115

In formulating the biodegradation, several simplifying assumptions were made towards

133

reducing model integration time and facilitating model calibration. Although ADM1 and ASM3

134

models consider multiple microbial populations, in the absence of Arctic-specific data on the

135

distribution of microbial populations, the proposed model only facultative heterotrophic (XB,H)

ACCEPTED MANUSCRIPT 7 and anaerobic (XB,AN) microorganisms. These populations are assumed to be present in all layers.

137

Also, only readily biodegradable soluble (Ss) and particulate (Xs) organic materials are

138

considered in the material balances, i.e., the current version of the model does not include

139

nitrogen and phosphorus balances.

RI PT

136

The aerobic and anaerobic liquid layers are assumed to have a saturation (maximal

141

attainable) concentration (Xmax) of suspended solids. At the end of each day of model integration,

142

the solids in excess of this saturation concentration are assumed to settle directly into the sludge

143

layer. Another simplifying assumption is that removal (or settling) of solids from the aerobic and

144

anaerobic layers to the sludge layer is instantaneous. Also, all dissolved components are assumed

145

to be uniformly distributed between the layers (instant mixing assumption). The instantaneous

146

settling and mixing assumptions can be justified by the fact that the lagoon performance is

147

simulated over a period of one year. Therefore the settling and mixing processes occurring at a

148

scale of several hours or even days may be considered instantaneous without any significant loss

149

of accuracy.

TE D

M AN U

SC

140

Material balances and kinetic equations describing the hydrolysis of particulate organic

151

materials and the growth of XB,H and XB,AN microbial populations on soluble organic materials

152

were adapted from ASM3 (Henze et al. 2000). In particular, multiplicative Monod kinetic

153

equations are used to describe both the aerobic and anaerobic (anoxic) growth of heterotrophic

154

biomass. For example, aerobic growth of the facultative heterotrophic biomass is described by

155

two multiplicative Monod kinetic terms representing available dissolved substrate (Ss) and

156

dissolved oxygen (SO) concentrations. The anaerobic heterotrophic biomass growth is also

157

described by two multiplicative terms representing growth on Ss and an oxygen inhibition

158

dependence, giving the following overall kinetic equation:

AC C

EP

150

ACCEPTED MANUSCRIPT 8 =

,



,

+

,

,

,

160

where

161

maximum specific growth rate of the heterotrophs,

162

constants, and

,

,

,

,

is the specific growth rate of the heterotrophic microorganisms,

,

,

,

and

is the inhibition constant.

,

,

,

,

,

,

,

) is described as

where

166

constant of anaerobic microorganisms and

,

is the maximum specific growth rate of the anaerobes, ,

M AN U

165

,

,

(1)

is the

SC

=

164

,

are the half saturation

The specific growth rate of anaerobic microorganisms (

163

167

,

RI PT

159

,

(2)

is the half-saturation

is the self-inhibition constant.

The rate of hydrolysis of particulate materials plays an important role in determining the

168

outcome of wastewater treatment, as it often represents the slowest biodegradation step. The

169

following expression, adapted from ASM3 model, is used to describe the specific hydrolysis rate

173 174 175

!" = #", $

(,

'%

TE D

172

hydrolysis:

,

')

,

,

&

,

*

+(

,

,

where #", is the maximum hydrolysis rate, #

EP

171

(!" ) in the proposed model, where both aerobic and anaerobic populations contribute to the

hydrolysis and - is the fraction of anaerobes(

AC C

170

+-

,

,

,

,

and #

),

,

(3) are the half saturation constants of

) contributing to hydrolysis.

It is to be noted that the kinetic parameters of the model are assumed to be temperature-

176

independent. This is because anaerobic microorganisms are known to adapt well to psychrophilic

177

conditions (Kashyap et al. 2003), therefore biodegradation by a complex consortium of

178

psychrophilic and mesophilic microorganisms can be described using such simplifying

179

assumption. Also, the rates of aerobic and anaerobic degradation of soluble organic materials are

ACCEPTED MANUSCRIPT 9 assumed to be proportional to the growth rates of the corresponding microbial populations.

181

Finally, the Arctic is considered to have near-arid conditions due to low precipitation rates.

182

Furthermore, all the lagoons examined for this study are man-made, which means that

183

surrounding embankments do not allow for runoff water to drain into the lagoon. Consequently,

184

dilution due to precipitation was not considered in the current version of the model. A detailed

185

description of these equations and the material balances of the model are provided in

186

Supplementary Information (Appendix B).

SC

RI PT

180

191 192 193 194 195 196

concentration in the lagoon (total - CODt, soluble - CODs and particulate - CODp) calculated as ./01 = ./07 =

2

34 2

34

∑, )

,

+

,

*6 ,

∑,)8 , + 8 , *6 ,

TE D

190

At each moment in time the performance of the WSP is described by the COD

(4) (5)

./09 = ./07 + ./01

(6)

where Vt is the combined volume of all layers (volume of the lagoon content), i is the

layer index, : = ;<, =., ;>, 8? corresponding to the aerobic, ice, anaerobic, and sludge layers.

EP

189

2.4 Numerical methods

First-order differential equations, corresponding to the model material balances, are

AC C

188

M AN U

187

197

numerically solved using ODE45 function of Matlab R2010 (Mathworks, Natick, MA, USA).

198

Model parameters are estimated using Matlab’s unconstrained optimization function fminsearch,

199

to minimize the root mean squared error (RMSE) between estimated and observed values of

200

total, soluble and particulate COD concentrations:

201

@A8
2

E

CI7 E ∑HM2 )FG,H − FG,HCK * & . L

(7)

ACCEPTED MANUSCRIPT 10

205 206 207 208 209 210 211 212

CK

represents the corresponding model output (estimated) values; k is the

index corresponding to COD components, k = N./09 , ./07 , ./01 O.

RI PT

204

observed) values; F

As will be described later, the sensitivity of model results to parameter values was used in

the calibration/validation process. The relative sensitivity PG of each model parameter Q with respect to each output variable FG was computed as: PG =

RSE RTU

∙S

TU

E

SC

203

Here, >G is the total number of measurements; F CI7 represents experimentally measured (or

(8)

where # = N./09 , ./07 , ./01 O and j is the parameter’s index.

M AN U

202

3. MODEL CALIBRATION AND VALIDATION

The proposed model was calibrated and validated using data collected at three arctic lagoons situated in the hamlets of Pond Inlet, Gjoa Haven, and Kugaruuk (Nunavut, Canada).

214

The three are engineered lagoons, with surface areas (when full) of about 40,000, 20,000, 13,000

215

m2, and average depths (when full) of about 1.6, 3 and 3 m, respectively. Available measurements at these lagoons included characterization of raw sewage and

EP

216

TE D

213

samples withdrawn from the lagoons (typically near the discharge point) during one year. Only

218

two to three sampling campaigns per year were carried out at each location, typically in early or

219

mid-summer as well as late-summer or early fall, towards decanting. Information provided

220

through these campaigns in include total and soluble COD and/or BOD concentration, and total

221

suspended solid (TSS) concentration. If total COD concentration values were unavailable, they

222

were estimated using BOD5 measurements by assuming a COD to BOD5 ratio of 2 (i.e., 1 mg L-1

223

of BOD5 is equivalent to 2 mg L-1 of COD). Biochemical measurements at Gjoa Haven and

224

Kugaruuk were provided by Environment and Climate Change Canada (ECCC, unpublished data

AC C

217

ACCEPTED MANUSCRIPT 11 on lagoon characterization from September 2009 to September 2010), while corresponding daily

226

temperatures, precipitation (rain and snow) and depth of snow on the ground were obtained from

227

the publicly available ECCC website and used for icecap formation/melting dates and thickness

228

computations. Biochemical data for Pond Inlet lagoon were taken from Ragush et al (2015).

229

These values represent averages of lagoon sampling between 2011 - 2014. Since the icecap

230

thickness and formation/melting dates in Pond Inlet were similar for all years, temperature values

231

retrieved from the ECCC website for 2009 – 2010 were used to calculate icecap properties.

SC

RI PT

225

232

234

3.1 Icecap formation

M AN U

233

The icecap formation and decay model described in Appendix A was used to estimate the beginning and end of ice formation as well as the icecap thickness. Temperature and recorded

236

snowfall information from existing Environment Canada weather stations close to lagoon

237

coordinates was used. Weather station identification is provided in Table 1. Figure 2 shows the

238

resulting ice thickness profiles calculated for the Pond Inlet, Gjoa Haven, and Kugaruuk lagoons

239

in Nunavut, Canada, based on the temperature profiles from September 2009 to September 2010.

240

In the fall of 2009 icecap started forming around September 15th to September 20th (relatively

241

soon after decant) in all three locations. Consequently, September 15th was considered as t = 0

242

(see Figure 2) for the model calibration calculations. The calculations of ice thickness were

243

carried out for 350 days. Water discharge was assumed to take place between days 350 – 365

244

(not shown in Figure 2) and effective average snow cover thickness was assumed to be 50% of

245

the “snow on the ground” values provided by Environment Canada data.

246 247

AC C

EP

TE D

235

As can be seen from Figure 2, at all three locations a near constant ice formation rate is estimated, resulting in a linear increase of ice thickness up to about day 200. Subsequently, a

ACCEPTED MANUSCRIPT 12 constant ice thickness is estimated between days 200-270, followed by a period of relatively fast

249

thaw that culminates with ice-free surface at about day 295. It should be noted that only 50% of

250

the accumulated snow on the ground (as recorded by Environment-Canada) was effectively

251

considered for heat-flux calculations, as a portion of the snow was assumed to have been blown

252

away due to winds.

253

RI PT

248

As mentioned earlier, the ice-formation model was not directly coupled with the lagoon model. Rather, these simulation results were used to determine ice formation parameters (start of

255

ice formation and thaw as well as constant ice formation and thaw rates), which were used as

256

inputs for the multi-layer lagoon model.

258 259

3.2 Biodegradation kinetics

M AN U

257

SC

254

Initially, model parameters related to heterotrophic activity were taken from ASM3 (Henze et al. 2000), while degradation rates under anaerobic conditions were assumed from (Lettinga et al.

261

2001). Also, initial value of the oxygen transfer coefficient (kL) was adapted from the work of

262

Chaturvedi et al. (Chaturvedi et al. 2014), assuming wind speeds between 1 and 8 m s-1 and an

263

average water temperature of 10 ºC. Table 2 provides the values of all model parameters and

264

Tables 3 and 4 describe influent wastewater composition and initial conditions used for model

265

integration.

EP

AC C

266

TE D

260

Considering the rather large number of model parameters and the limited number of available

267

measurements (only 2-3 sets of COD and TSS measurements were available for each lagoon,

268

with no biomass measurements at all), it was deemed best to identify a small subset of model

269

parameters with the largest influence on model outputs and subsequently calibrate the model

270

using only this limited subset, while all other parameters are kept constant at values obtained

ACCEPTED MANUSCRIPT 13 271

from the literature. The most influential parameters were identified through a sensitivity analysis,

272

using Eq. 8. Table 4 summarizes results of this sensitivity analysis. Based on these results, a

273

subset of three parameters was selected for estimation, namely the maximum specific growth rate

274

of anaerobic microorganisms (

275

oxygen transfer coefficient (kL).

), the rate of solid organic matter hydrolysis (#" ), and the

RI PT

276

,

Parameter estimation (model calibration) was performed using data of two lagoons, Pond , W,

#" , and kL that minimized RMSE

Inlet and Gjoa Haven, whereby values of parameters

278

(Eq. 7) were identified (all other parameters were held constant). Figure 3 provides the resulting

279

profiles of total, soluble, and particulate CODs for the Pond Inlet and Gjoa Haven lagoons with

280

the calibrated parameters. Table 6 summarizes the estimated values of these parameters. The

281

model was subsequently validated by applying the estimated parameter values to the Kugaruuk

282

lagoon. Figure 4 illustrates a comparison between estimated and observed values of COD and

283

TSS measurements in the Kugaaruk lagoon.

M AN U

TE D

284

SC

277

It can be seen that in all cases the model predicted a COD accumulation during winter months followed by a relatively fast COD degradation during the ice-free period. While the

286

predicted trends qualitatively agreed with the observed values, the RMSE value for Kugaruuk

287

validation was higher (as expected) than for the lagoons used for parameter estimation (Table 5).

288

It was nonetheless deemed that the model showed a reasonable capability to correctly predict the

289

measured trends in total and dissolved COD concentrations.

AC C

290

EP

285

As mentioned earlier, the small number of measured data available for calibration

291

presented a significant challenge to the confidence with which the results could be evaluated.

292

Consequently, an additional calibration method was undertaken whereby the estimation

293

procedure was repeated for each lagoon individually. In this case however, the oxygen transfer

ACCEPTED MANUSCRIPT 14 , W

and #", values were re-

294

coefficient was kept at the estimated value of 6.4 m d-1 , while

295

estimated for each lagoon. As shown in Table 6, this approach resulted in smaller RMSE values

296

and also demonstrated the differences between the lagoons in terms of

297

Since parameter sensitivity analysis showed high impact of these parameters on the model

298

outputs, the differences between the lagoons are likely represented by the differences between

299

lagoon-specific values.

RI PT

and #", values.

It is acknowledged that the scarcity of field data, especially in the winter, limits the

SC

300

, W

number of estimated model parameters and the accuracy of the parameter estimation procedure.

302

The authors are currently involved in an effort to augment available data through year-round

303

sampling of a northern lagoon near Yellowknife (NWT, Canada). While it is too early to report

304

additional results in this study, future work will hopefully be based on a much larger dataset.

M AN U

301

306 307

4. DISCUSSION

TE D

305

While the dearth of calibration data limits the confidence in the accuracy of the results, some important, and not altogether expected, trends on the performance of facultative sewage

309

lagoons in the Arctic can be discerned from the model outputs. For example, a clear trend

310

emerges, whereby soluble COD concentration increases consistently during winter months

311

(Figures 3, 4). This can be explained by the fact that during winter, particulate COD settles and

312

continuously undergoes hydrolysis in the sludge layer of the lagoon, albeit at a low rate, while at

313

the same time a relatively slow anaerobic degradation acts to reduce soluble COD concentration.

314

When the icecap melts away, fast reduction in soluble COD concentration can be observed, as

315

aerobic activity initiates and increases in the top layer. The aerobic growth of suspended

316

heterotrophic biomass is accompanied by biomass settling, since the suspended biomass

AC C

EP

308

ACCEPTED MANUSCRIPT 15 concentration in the aerobic layer is limited by the maximum attainable density (Xmax). The

318

aerobic biomass in the sludge layer contributes to hydrolysis of solids according to Eq. (3).

319

Notably, the hydrolysis rate is not limited by the absence of oxygen in the sludge. As a result,

320

the overall rate of hydrolysis accelerates during summer.

321

RI PT

317

Another important observation can be made that although the degradation of organic materials under anaerobic conditions is slow (compared to aerobic degradation), the anaerobic

323

activity provides a significant contribution to the overall COD removal in the lagoon due to the

324

fact that this activity takes place throughout the year and also due to the high sludge density. This

325

can be further demonstrated by comparing different scenarios of COD degradation in the shallow

326

Pond Inlet lagoon (Figure 5A) and the deep Kugaruuk lagoon (Figure 5B). In one scenario the

327

COD degradation is simulated with the anaerobic rate of COD degradation set to zero (only

328

aerobes), while in the other scenario the aerobic degradation rate is set to zero (only anaerobes).

329

It is clear that the absence of anaerobic biodegradation throughout the winter results in COD

330

accumulation in both lagoons. If zero anaerobic activity is assumed, the rate of aerobic

331

degradation required to achieve observed COD concentrations at the end of summer would have

332

to be increased to an unrealistically high value of approximately 10 d-1, which is not likely given

333

that water temperature in the lagoons does not typically exceed 10-12oC during the summer.

334

Clearly, only a combination of both aerobic and anaerobic degradation activities results in the

335

COD profiles matching the field measurements.

M AN U

TE D

EP

AC C

336

SC

322

Interestingly, the model predicts different COD profiles during winter for each lagoon.

337

Pond Inlet lagoon is relatively shallow (about 1.6 m average depth when full), in which the ice

338

thickness can reach a maximum of about 1.4 m (Figure 2). This results in an expected gradual

339

increase of total and soluble COD concentrations throughout the winter due to the absence of

ACCEPTED MANUSCRIPT 16 biological activity in the ice layer. Although, anaerobic degradation continues in the sludge layer,

341

it is not sufficient to considerably reduce soluble CODs (Figure 3). This prediction agrees well

342

with the high COD value measured at the beginning of summer, right after the complete melting

343

of the ice cover. This accumulation of CODs is followed by a relatively fast COD reduction in

344

the summer. It can be seen that even in the shallow Pond Inlet lagoon anaerobic degradation

345

contributes to the overall COD removal as the comparison of different scenarios shown in Figure

346

5A indicates.

SC

RI PT

340

The predicted COD profiles in the much deeper Gjoa Haven, and Kugaruuk lagoons followed

348

a somewhat different pattern. In these lagoons a fast increase in COD concentrations is predicted

349

during the first 50-60 days, followed by slow COD concentration decrease throughout the winter,

350

due to anaerobic activity in the unfrozen anaerobic and sludge layers under the ice. The rate of

351

COD degradation is predicted to accelerate during the summer (Figures 3 and 4) due to aerobic

352

activity. It is noted (again) that only COD trends in the summer months can be confirmed by

353

observed data. Profiles of particulate CODs in all lagoons were similar with a somewhat higher

354

concentration during summer. Fast growth of aerobic bacteria in the upper oxygenated layer

355

during summer months results is accompanied by biomass settling, which contributes to

356

particulate CODs as well as to sludge accumulation, although it increases the rate of hydrolysis

357

according to model equation (3). Notably, solids hydrolysis does not require dissolved oxygen

358

and therefore proceeds in the sludge layer, where the biomass concentration is the highest.

TE D

EP

AC C

359

M AN U

347

A detailed look at the model outputs and the parameter estimation results suggests that

360

the difference in COD profiles can be explained by the differences in lagoon depths. The shallow

361

Pond Inlet lagoon is frozen during 260 days of the year (most of this period frozen through,

362

except the sludge layer), which significantly limits anaerobic activity during winter months.

ACCEPTED MANUSCRIPT 17 During the summer months a high degradation rate is achieved due to higher ratio of aerobic to

364

anaerobic activity. In Pond Inlet, an aerobic layer thickness of 40 cm corresponds to

365

approximately 25% of this lagoon depth, while in two other lagoons the aerobic layer comprises

366

less than 15% of the total lagoon depth. Consequently, aerobic degradation provides a relatively

367

higher contribution to COD removal at the Pond Inlet lagoon compared to two other, deeper

368

lagoons (Figure 5).

RI PT

363

Gradual degradation of COD throughout winter in the Gjoa Haven and Kugaruuk lagoons

370

is attributed to the existence of a liquid anaerobic layer throughout the winter months (Figures 3,

371

4), which supports anaerobic activity. As a result, a moderate COD concentration decrease is

372

predicted, despite the presence of the icecap.

M AN U

373

SC

369

The importance of considering anaerobic microbial activity in northern lagoons appears to be an immediate conclusion, with important implications for future lagoon design. It can be

375

surmised that some combination of deep and shallow lagoons is likely to provide the best results

376

in terms of degradation of organics in the sewage. In fact, there is some evidence that multi-cell

377

lagoons in the Arctic are functioning well. It appears that a deep first lagoon is likely to provide

378

substantial anaerobic degradation in winter, while a shallow and wide second lagoon will provide

379

substantial aerobic degradation in summer.

EP

Another important consideration in evaluating Arctic lagoon performance could be

AC C

380

TE D

374

381

related to the contribution of algae blooms to COD removal. This attribution is supported by

382

Ragush et al (2018) and Schmidt et al (2017), who examined oxygen dynamics and organic

383

degradation for lagoons during the ice-free period and demonstrated algae’s contribution in

384

supporting aerobic degradation. While the current version of the model does not consider algae

385

growth, it appears to be indirectly accounted for by the oxygen transfer coefficient (kL). During

ACCEPTED MANUSCRIPT 18 long summer days algae are expected to supply oxygen to the top (aerobic) layer of the lagoon.

387

Indeed, parameter estimation yielded a kL value of 6.4, which is somewhat higher than a range of

388

values suggested in literature (0.9 – 6.0 m d-1), as indicated in Table 2. In addition to increased

389

oxygen supply, both phototrophic and heterotrophic growth of algae are likely to contribute to

390

solids accumulation, especially after the end of summer, when temperatures become lower and

391

days become shorter, before the formation of icecap. Finally, consumption of dissolved methane

392

by aerobic methanotrophic bacteria could be another factor affecting both dissolved oxygen and

393

methane concentrations throughout summer. Future versions of the Arctic lagoon model need to

394

be extended to account for these phenomena as well as for nitrogen and phosphorus removal at

395

low temperatures.

M AN U

SC

RI PT

386

396

398

5. CONCLUSION

This study describes a mathematical model capable of simulating the biodegradation of

TE D

397

wastewater in arctic WSPs. The model accounts for several unique features of Arctic lagoons,

400

including the periodic nature of lagoon operation and the existence of icecap for most of the

401

lagoon operation cycle. The proposed model is calibrated and validated using limited observed

402

data, obtained at three arctic lagoons. Within the constraints of limited calibration data and

403

previously noted limitations of the parameter estimation procedure, the differences between

404

predicted and measured COD and TSS concentration profiles (in terms of mean squared error)

405

are deemed to be indicative of an acceptable qualitative agreement between them. Analysis of

406

predicted COD concentration profiles suggests significant anaerobic biodegradation of organic

407

materials throughout winter months. The next phase of model development will incorporate field

408

data from a northern sewage lagoon with wintertime measurements of BOD, COD, TSS and

AC C

EP

399

ACCEPTED MANUSCRIPT 19 temperatures at different depths beneath the ice. Investigation into the contribution of algae and

410

enhancement of the ice thickness model are also planned. Additional data are expected to nuance

411

and confirm key assumptions as well as improve the model as a tool for analyzing and improving

412

Arctic lagoon design.

RI PT

409

413

Acknowledgement

415

The authors wish to acknowledge Environment and Climate Change Canada (ECCC) and in

416

particular Ms. Shirley Anne Smyth, who provided unpublished data that were invaluable for the

417

initial calibration of the model. Also Mr. Babesh Roy of the Government of Nunavut, who

418

facilitated data sharing on the Pond Inlet lagoon.

AC C

EP

TE D

M AN U

SC

414

ACCEPTED MANUSCRIPT 20 Notations

;

Specific interfacial area, m-1

BOD

Y

Biological oxygen demand, mg L-1

Y

Decay rate of Decay rate of

COD

Chemical oxygen demand, mg L-1

CSTR

Continuous stirred tank reactor

ℎ

Fraction of

Z[

,

*, d-1 *, d-1

yielding

,%

Lagoon’s depth, m

, ,

Self-inhibition coefficient for

,

, mg L-1

M AN U

,

Oxygen mass transfer coefficient**, m d-1

Half-saturation coefficient of 8 for Inhibition coefficient of 8 for

,

,

Half-saturation coefficient of 8 for

,

Half-saturation coefficient of

,

,

Half-saturation coefficient of 8 for

for

Maximum specific hydrolysis rate*, d-1

, ,

,

MSE

^,W

Mean squared error

^9_

Input flow rate, L d-1

8

Concentration of soluble oxygen, mg L-1

8

87

,

cd b,W,

cd beWK

6

WSP

, mg L-1

, mg L-1

, mg L-1

Transfer flow rate from layer : to layer a, L d-1

EP

8

,→H

, mg L-1

, mg L-1

TE D

]

#"

,

SC

X

RI PT

Lagoon’s surface, m2

Maximum value for 8 , mg L-1

AC C

419

Concentration of soluble inert materials, mg L-1 Concentration of biodegradable soluble substrate, mg L-1

Start of ice formation, d End of ice thaw, d Lagoon’s volume, L

Wastewater stabilization pond

ACCEPTED MANUSCRIPT 21 ,

Concentration of heterotrophic bacteria, mg L-1

,



d,

Concentration of anaerobic bacteria, mg L-1 Concentration of particulate inerts, mg L-1

Solids saturation density in sludge, mg L-1

c

Solids saturation density in ice, mg L-1

f,

Yield factor of 8 for

g

layer height, m

g

h

d

-"

Solids settling, %

c

Period of ice formation/thaw, d

Correction factor for hydrolysis by Correction factor anoxic growth , ,

Maximum specific growth rate of Maximum specific growth rate of

;>

klm kno p ℎ

: or a

Aerobic layer

Anaerobic layer

AC C

;<

EP

Indices

=.

, mg/mg

mg/mg

Maximum thickness of ice layer, m

-j

=

,

,

Maximum thickness aerobic layer, m

c

Δb

mg/mg

,

TE D

g

Yield factor of 8 for

,

,

M AN U

Yield factor of 8 for

SC

Concentration of biodegradable particulate substrate, mg L-1

f, f

RI PT

Solids saturation density in liquid, mg L-1

Final

Experimental values

Heterotrophic microbial population

Related to hydrolysis Inert Ice layer Any of the layers in the model

, ,

,

*

*, d-1 *, d-1

ACCEPTED MANUSCRIPT 22 :l:

qXn

Initial

8?

P:q

Sludge layer or related to the soluble substrate

Maximum

/

Oxygen, 8 , mg L-1

b!

Simulated values Transfer between layers

SC

420

Influent

RI PT

:l

421

AC C

EP

TE D

M AN U

422

ACCEPTED MANUSCRIPT 23 REFERENCES

424

Aldana, G.J., Lloyd, B.J., Guganesharajah, K. and Bracho, N. (2005) The development and

425

calibration of a physical model to assist in optimising the hydraulic performance and design of

426

maturation ponds. Water Sci. Technol. 51(12), 173-181.

427

Ashton, G.D. (1983) Predicting Lake Ice Decay, Cold regions research and engineering lab

428

Hanover NH.

429

Ashton, G.D. (1986) River and lake ice engineering, Water Resources Publication.

430

Ashton, G.D. (1989) Thin ice growth. Water Resour. Res. 25(3), 564-566.

431

Banks, C.J., Koloskov, G.B., Lock, A.C. and Heaven, S. (2003) A computer simulation of the

432

oxygen balance in a cold climate winter storage WSP during the critical spring warm-up period.

433

Water Sci. Technol. 48(2), 189-196.

434

Batstone, D.J., Keller, J., Angelidaki, I., Kalyuzhnyi, S.V., Pavlostathis, S.G., Rozzi, A.,

435

Sanders, W.T.M., Siegrist, H. and Vavilin, V.A. (2002) The IWA anaerobic digestion model no

436

1 (ADM1). Water Sci. Technol. 45(10), 65-73.

437

Beran, B. and Kargi, F. (2005) A dynamic mathematical model for wastewater stabilization

438

ponds. Ecol. Model. 181(1), 39-57.

439

Borja, R., González, E., Raposo, F., Millán, F. and Martín, A. (2002) Kinetic analysis of the

440

psychrophilic anaerobic digestion of wastewater derived from the production of proteins from

441

extracted sunflower flour. J. Agric. Food Chem. 50(16), 4628-4633.

442

Chaturvedi, M.K.M., Langote, S.D., Kumar, D. and Asolekar, S.R. (2014) Significance and

443

estimation of oxygen mass transfer coefficient in simulated waste stabilization pond. Ecol. Eng.

444

73, 331-334.

445

Chen, Y.R. and Hashimoto, A.G. (1978) Kinetics of methane fermentation, pp. 269-282.

AC C

EP

TE D

M AN U

SC

RI PT

423

ACCEPTED MANUSCRIPT 24 Conrad, R. (1999) Contribution of hydrogen to methane production and control of hydrogen

447

concentrations in methanogenic soils and sediments. FEMS Microbiol Ecol 28, 193-202.

448

Dochain, D., Gregoire, S., Pauss, A. and Schaegger, M. (2003) Dynamical modelling of a waste

449

stabilisation pond. Bioprocess. Biosyst. Eng. 26(1), 19-26.

450

Gu, R. and Stefan, H.G. (1995) Stratification dynamics in wastewater stabilization ponds. Water

451

Res. 29(8), 1909-1923.

452

Henze, M., Gujer, W., Mino, T. and van Loosdrecht, M.C.M. (2000) Activated sludge models

453

ASM1, ASM2, ASM2d and ASM3, IWA publishing.

454

Houweling, C.D., Kharoune, L., Escalas, A. and Comeau, Y. (2005) Modeling ammonia removal

455

in aerated facultative lagoons. Water Sci. Technol. 51(12), 139-142.

456

Hunik, J.H., Bos, C.G., Hoogen, M.P., DeGooijer, C.D. and Tramper, J. (1994) Co-immobilized

457

Nitrosomonas europea and Nitrobacter agilis cells: validation of a dynamic model for

458

simultaneous substrate conversion and growth in k-carrageenan gel beads. Biotechnology and

459

Bioengineering 43, 1153-1163.

460

Johnson, K., Prosko, G. and Lycon, D. (2014) The challenges with mechanical wastewater

461

systems in the far north, Regina, SK.

462

Jupsin, H. and Vasel, J.L. (2007) Modelisation of the contribution of sediments in the treatment

463

process case of aerated lagoons. Water Sci. Technol. 55(11), 21-27.

464

Kashyap, D.R., Dadhich, K.S. and Sharma, S.K. (2003) Biomethanation under psychrophilic

465

conditions: a review. Biores. Technol. 87, 147-153.

466

Kayombo, S., Mbwette, T.S.A., Mayo, A.W., Katima, J.H.Y. and Jorgensen, S.E. (2000)

467

Modelling diurnal variation of dissolved oxygen in waste stabilization ponds. Ecol. Model.

468

127(1), 21-31.

AC C

EP

TE D

M AN U

SC

RI PT

446

ACCEPTED MANUSCRIPT 25 Langergraber, G., Rousseau, D.P.L., García, J. and Mena, J. (2009) CWM1: a general model to

470

describe biokinetic processes in subsurface flow constructed wetlands. Water Sci. Technol.

471

59(9), 1687-1697.

472

Lettinga, G., Rebac, S. and Zeeman, G. (2001) Challenge of psychrophilic anaerobic wastewater

473

treatment. Trends Biotechnol. 19(9), 363-370.

474

Martínez, F.C., Cansino, A.T., García, M.A.A., Kalashnikov, V. and Rojas, R.L. (2014)

475

Mathematical analysis for the optimization of a design in a facultative pond: indicator organism

476

and organic matter. Math. Probl. Eng. 2014.

477

Peng, J.F., Wang, B.Z., Song, Y.H. and Yuan, P. (2007) Modeling N transformation and removal

478

in a duckweed pond: model development and calibration. Ecol. Model. 206(1), 147-152.

479

Perry, R. and Green, D. (1984) Perry's Chemical Engineers' Handbook, McGraw-HiII, Inc.

480

NewYork.

481

Ragush, C.M., Gentleman, W.C., Truelstrup-Hansen, L. and Jamieson, R.C. (2018) Operational

482

Limitations of Arctic Waste Stabilization Ponds: Insights from Modeling Oxygen Dynamics and

483

Carbon Removal. Environ Eng 144, 12.

484

Ragush, C.M., Schmidt, J.J., Krkosek, W.H., Gagnon, G.A., Truelstrup-Hansen, L. and

485

Jamieson, R.C. (2015) Performance of municipal waste stabilization ponds in the Canadian

486

Arctic Ecological Eng. 83, 413-421.

487

Rauch, W., Vanhooken, H. and Vanrolleghem, P.A. (1999) A simplified mixed-culture biofilm

488

model. Water Research 33, 2148-2162.

489

Sah, L., Rousseau, D.P.L. and Hooijmans, C.M. (2012) Numerical modelling of waste

490

stabilization ponds: where do we stand? Water Air Soil Pollut. 223(6), 3155-3171.

AC C

EP

TE D

M AN U

SC

RI PT

469

ACCEPTED MANUSCRIPT 26 Sah, L., Rousseau, D.P.L., Hooijmans, C.M. and Lens, P.N.L. (2011) 3D model for a secondary

492

facultative pond. Ecol. Model. 222(9), 1592-1603.

493

Salter, H.E., Ta, C.T., Ouki, S.K. and Williams, S.C. (2000) Three-dimensional computational

494

fluid dynamic modelling of a facultative lagoon. Water Sci. Technol. 42(10-11), 335-342.

495

Schmidt, J.J., Gagnon, G.A. and Jamieson, R.C. (2017) Predicting microalgae growth and

496

phosphorus removal in cold region waste stabilization ponds using a stochastic modeling

497

approach. Environ Sci Pollut Res, 12.

498

Stephenson, R.J., Patoine, A. and Guiot, S.R. (1999) Effects of oxygenation and upflow liquid

499

velocity on a coupled anaerobic/aerobic reactor system. Water Research 33, 2855-2863.

500

Sweeney, D.G., Cromar, N.J., Nixon, J.B., Ta, C.T. and Fallowfield, H.J. (2003) The spatial

501

significance of water quality indicators in waste stabilization ponds-limitations of residence time

502

distribution analysis in predicting treatment efficiency. Water Sci. Technol. 48(2), 211-218.

503

Sweeney, D.G., Nixon, J.B., Cromar, N.J. and Fallowfield, H.J. (2005) Profiling and modelling

504

of thermal changes in a large waste stabilisation pond. Water Sci. Technol. 51(12), 163-172.

505

Tartakovsky, B. and Guiot, S.R. (1997) Modeling and analysis of layered stationary anaerobic

506

granular biofilms. Biotechnology and Bioengineering 54, 122-130.

507

Tartakovsky, B., Manuel, M.F. and Guiot, S.R. (2006) Degradation of trichloroethylene in a

508

coupled anaerobic-aerobic bioreactor : modeling and experiment. Biochemical Engineering

509

Journal 26(1), 72-81.

510

Toprak, H. (1994) Empirical modelling of sedimentation which occurs in anaerobic waste

511

stabilization ponds using a lab‐scale semi‐continuous reactor. Environ. Technol. 15(2), 125-134.

512

van der Berg, L. (1977) Effect of temperature on growth and activity of a methanogenic culture

513

utilising acetate. Can. J. Microbiol. 23(7), 898-902.

AC C

EP

TE D

M AN U

SC

RI PT

491

ACCEPTED MANUSCRIPT 27 Yates, C.N., Wootton, B.C. and Murphy, S.D. (2012) Performance assessment of arctic tundra

515

municipal wastewater treatment wetlands through an arctic summer. Ecological Eng. 44, 160-

516

173.

RI PT

514

517 518

SC

519 520

AC C

EP

TE D

M AN U

521

ACCEPTED MANUSCRIPT 28 Table 1. Climate station identifying information. Station

Coordinates

Elevation (m)

Latitude

Longitude

Gjoa Haven, NU

68.64

-95.85

42.3

Kugaaruk, NU

68.54

-89.8

15.5

Pond Inlet, NU

72.69

-77.97

61.6

523

Climate

WMO

TC

2302340

71363

ZHK

2303092

2403202

AC C

EP

TE D

M AN U

SC

524

Identifiers

RI PT

522

71095

YBB YIO

ACCEPTED MANUSCRIPT 29 Table 2. Values of the lagoon model parameters. Parameter Value Notes Source 1 3 Assumed (Henze et al. 2000) , 1 0.09 – 0.12 Calibrated (Batstone et al. 2002) , (van der Berg 1977) (Chen and Hashimoto 1978) (Borja et al. 2002) 0.051 Assumed (Henze et al. 2000) Y 0.021 Assumed (Batstone et al. 2002) Y 20 Assumed (Henze et al. 2000) , 28 Assumed (Batstone et al. 2002) , (Langergraber et al. 2009) (Borja et al. 2002) 0.2 Assumed (Henze et al. 2000) , (Banks et al. 2003) 0.01 Assumed (Henze et al. 2000) , `1 0.10 - 0.22 Calibrated (Henze et al. 2000) #" 0.11 Assumed (Henze et al. 2000) -" (Langergraber et al. 2009) 0.2 Assumed (Henze et al. 2000) , 200 Assumed , 12 Assumed (Perry and Green 1984) 8 , (Sah et al. 2011) 2 6.4 Calibrated (Chaturvedi et al. 2014) ] 0.32 0.6 Measured X 0 Assumed (Henze et al. 2000) -j 1.58 Assumed (Henze et al. 2000) f, 1.3 Assumed (Batstone et al. 2002) f, 1.72 Assumed (Henze et al. 2000) f, 0.08 Assumed (Henze et al. 2000) Z[ 11,093Measured this work ; 36,686 1.65-3.13 Measured this work ℎ c 10 Calculated this work b,W, c 295 Calculated this work beWK c 30 Calculated this work Δb c 1.2 - 1.5 Calculated this work g 90 Assumed this work h d, 32.0-70.3 Measured this work ] 8,000 Assumed this work c infinite Assumed this work d 0.4 Assumed this work g 1 Set to zero in the ice layer 2 Only considered for the aerobic layer

Range 0.6 – 13.2 0.05 – 8

RI PT

525

AC C

EP

TE D

M AN U

SC

0.05-1.6 5-225 28 - 150

526 527

0.01 – 0.2

0.01 – 0.03 1-3 0.1-0.4

0.01 – 0.2 7-15 0.86 – 6.04 0.6 –.0 1.33 – 2.63 0.08 -

ACCEPTED MANUSCRIPT 30 Table 3. Influent wastewater composition.

Heterotrophic biomass ,

Anaerobic biomass Soluble biodegradable substrate Particulate biodegradable substrate Dissolved oxygen Particulate inert substrate Soluble inert substrate

87,,W 8

,,W

,,W

,,W

8 ,,W

Value 5

mg L-1

5

mg L-1

10101; 7892; 7753

Measured

mg L-1

4601; 3912; 3663

Measured

mg L-1 mg L-1 mg L-1

Assumed Assumed

0

Assumed

0

Assumed

0

Assumed

138.11; 128.22;

3

m d

TE D EP

Notes

mg L-1

-1

85.13

Pond Inlet lagoon, 2Gjoa Haven lagoon, and 3Kugaruuk lagoon

AC C

530

1

,W

^,W

Flow rate

529

, ,W

Dimension

RI PT

Notation

SC

WW component

M AN U

528

Measured

ACCEPTED MANUSCRIPT 31 531

Table 4. Initial conditions used for lagoon model integration. Layer Dimension

Notes

87



mg L-1

10

0

0

mg L-1

0

10

100

1

mg L

260 ; 200

mg L-1

3.5

mg L-1

0

mg L-1

0

m

g

d

1

2,3

260 ; 200

641; 502,3

1

2,3

260 ; 200 ]

0

Assumed

0

Assumed

0

Assumed2,3

Measured1

0 Assumed2,3

0

0

0

Measured

0

0

0

Assumed

0

0

0

Assumed

0

0.02

0

Assumed

Pond Inlet lagoon, 2Gjoa Haven lagoon, and 3Kugaruuk lagoon

EP

1

641; 502,3

AC C

532 533 534 535

g

2,3

TE D

8

Ice

Measured1

-1

mg L-1 8

Sludge

M AN U

,

,

Anaerobic

RI PT

Aerobic

SC

Parameter

ACCEPTED MANUSCRIPT 32 Table 5. Parameter sensitivity values with respect to CODt, CODp, and CODs variations.

, ,

,

,

-6.8 -256.9 5.5 0.9 13.9 0.2 -0.1 285.4 14.5 -0.1 -997.9 -3.8

TE D EP

g

,

] d

-12.6 710.4 48.8 57.8 -2.6 0.8 0.1 -5484.9 85.1 -0.4 -430.5 1.0

AC C

538 539

#" -"

stuvx

stuvs

-15.5 968.8 59.9 69.8 -5.4 1.0 0.2 -6760.8 124.0 -0.5 -387.3 2.3

RI PT

Y Y

,

stuvw

SC

Parameter

M AN U

536 537

ACCEPTED MANUSCRIPT 33 540

Table 6. Model calibration and validation results. Notations: PI – Pond Inlet lagoon, GH - Gjoa

542

Haven, and KU – Kugaruuk.

543

PI and GH

kh

RMSE

Notes

6.40 0.043 0.07

197

calibration

KU

6.40 0.043 0.07

195

validation

PI

6.40 0.100 0.11

63

re-calibration

GH

6.40 0.220 0.09

69

re-calibration

KU

6.40 0.22

102

re-calibration

µmax,AN

0.12

544 545 546 547 548

AC C

EP

TE D

549 550

SC

KL

M AN U

Lagoon

RI PT

541

ACCEPTED MANUSCRIPT 34 551 552

^,W ^,W]

]

(b)

g

Anaerobic

Sludge

g ] (b)

554 555 556

558 559

565

^,W

^,W ^,W]

^9_c↔

Ice

g c (b)

Anaerobic

]

Sludge

g

(b)

g ] (b)

AC C

561

564

=y

EP

560

563

d

(B) Cold period

TE D

557

562

g

^,Wc

M AN U

553

^9_d→

Aerobic

/L

RI PT

^,W

(A) Warm period

SC

^,Wd

566

Figure 1. Schematic diagram of the layers considered in the model during (A) warm and (B)

567

cold periods.

ACCEPTED MANUSCRIPT 35 568 569

Formation

Saturation Thaw No ice

RI PT

Pond Inlet Gjoa Haven Kugaaruk

1.2

0.8

SC

Ice thickness (m)

1.6

0.4

0.0 0

50

571 572

578 579

250

300

350

EP

577

200

AC C

576

150

TE D

573

575

100

295

Time (days)

570

574

270

M AN U

200

580

Figure 2. Icecap profile simulated using temperature data collected in 2010. Ice thickness values

581

are presented from September 15th 2009 (day 0) until August 31st ,2010 (day 350).

582 583

ACCEPTED MANUSCRIPT 36 584 585 C Gjoa Haven

A Pond Inlet

586

RI PT

587 588 589

SC

590 B

D

591

M AN U

592 593 594 595

599 600 601 602

EP

598

AC C

597

TE D

596

603

Figure 3. Model calibration results showing a comparison of model outputs with experimental

604

data and predicted layer thickness for (A, B) Pond Inlet and (C, D) Gjoa Haven lagoons.

605

ACCEPTED MANUSCRIPT 37 606 607 A Kugaaruk

608

RI PT

609 610 611

613

B

M AN U

614 615 616 617

622 623

EP

621

AC C

620

TE D

618 619

SC

612

624

Figure 4. Model validation results showing a comparison of model outputs with experimental

625

data (A) and predicted layer thickness (B) for Kugaruuk lagoon.

626 627 628

ACCEPTED MANUSCRIPT 38 629 630 631 A Pond Inlet

B Kugaaruk

RI PT

632 633 634

SC

635 636

M AN U

637 638 639 640

644 645 646 647

EP

643

AC C

642

TE D

641

648

Figure 5. A comparison of COD concentration profiles calculated with and without contribution

649

of aerobic and anaerobic degradation for the Pond Inlet (A) and Kugaaruk (B) lagoons.

650

ACCEPTED MANUSCRIPT 39 651

SUPPLEMENTARY INFORMATION

652

APPENDIX A. Modelling icecap formation and decay

653

The following equation was used to describe ice formation. This equation is based on work by Ashton (Ashton 1986, 1989) with modifications to account for the heat flux from the

655

relatively warm sewage inflow to the ice layer and for the insulating properties of snow cover: ∆b ~ −~ } − p€, (~€ − ~ )• ℎ ℎ 1 |? % , + 7+p & #, #7

(A1)

SC

∆ℎ =

RI PT

654

where

657

∆h = incremental growth of ice layer thickness [m] over time-step ∆t

658

ki = thermal conductivity of ice [W m-1 oC-1] (taken as 2.3 (corresponding to temp. of

659

M AN U

656

about -10 oC) for this research)

ks = thermal conductivity of snow [W m-1 oC-1] (taken as 0.4 for this research)

661

Tm = temperature at the ice water interface (assumed 0oC)

662

Tw = temperature of the water under the ice

663

Ta = temperature of air above the ice [oC]

664

Ha = thermal resistance (or heat transfer coefficient) between the uppermost surface and

666 667 668 669 670 671

EP

the air above it (i.e., when hs = 0 then Ha = Hia; and when hs > 0 then Ha = Has)

AC C

665

TE D

660

Hia = bulk heat transfer coefficient between the top surface temperature of the ice and the

air temperature above the ice (taken as 25 [W m-2 oC-1] for this research

Hsa = bulk heat transfer coefficient between the top surface temperature of snow and the

air temperature above the snow (taken as 15 [W m-2 oC-1] for this research) Hwi = bulk heat transfer coefficient between the water and the ice (taken as 2.19 [W m-2 o -1

C ] for this research

ACCEPTED MANUSCRIPT 40 672 673

ρ=

ice density [kg/m-3] (taken as 919 (corresponding to temp. of about -10 oC) for this research)

L = heat of fusion [J/kg-1] (taken as 3.3355•105 for this research)

675

hi = thickness of the ice cover at time-step ∆t

676

hs = thickness of the snow cover at time-step ∆t

677

When ∆t is taken as one day (or 86,400 seconds) the ice layer thickness becomes a function of freezing degree-days.

SC

678

RI PT

674

The following equation was used to model the decay of ice over a sewage lagoon. This

680

equation is also based on work by Ashton (Ashton 1983) with modifications to account for the

681

heat flux from the relatively warm sewage inflow to the ice layer: ∆ℎ =

∆b ‚−p, (~ − ~ ) − p€, (~€ − ~ )ƒ |?

(A2)

Here too, when ∆t is taken as one day the ice layer thickness becomes a function of

TE D

682

M AN U

679

683

melting degree-days. In the absence of real field data, the water temperature Tw is assumed to

684

have a sinusoidal form over the year with an amplitude of ATw and a period of 365 days. The

685

average temperature in day d is therefore estimated as:

EP

~„… = ~„† ‡ + 0.5; B‹ Œ1 + •ŽP

687 688 689

(A3)

where d =1 is the day at which decanting ends. Minimum water temperature ~„† ‡ was assumed

AC C

686

2•m “ (m = 1, 2, … 365) 365

to be 3oC. The amplitude can be expected to be in the range of 4 to 12 degrees.

ACCEPTED MANUSCRIPT 41 690

APPENDIX B. Biodegradation model

691

The following material balance equations are used to describe sewage biodegradation. The

m

,

mb

,

=

8, 8, 8, , , $ + $ + $ +$ + ,, + , j , , , , , , ™ššššššššššš›šššššššššššœ ™ššššššššššš›šššššššššššœ , +8 , +8 , +8 , +8 Aerobic growth of heterotrophic biomass

−

Y ™š›šœ ,

Decay of heterotrophic biomass

2. Anaerobic biomass (•˜–,•ž ) ,

mb

=

Flow in

Transport term

8, , , , $ +$ +$ + ,, − Y , ™š š›š šœ , , , , + 8 + 8 + , , , ™ššššššššššššššššš›šššššššššššššššššœ , Decay of Growth of anaerobic biomas

, ^,W ^ H + , ) , ,,W − , , * + 9_, ) , − , , * ™ššš ššš›ššššššœ ™ššššš›šššššœ 6 6 H→,

TE D

,

, ^,W ^ H + , ) , ,W − , , * + 9_, ) , − , , * ™š 6 šššš›šššššœ ™š 6 šššš›šššššœ

Flow in

EP

m

Anoxic growth of heterotrophic biomass H→,

,

AC C

694

RI PT

1. Heterotrophic biomass (•˜–,— )

SC

693

balances are composed for each layer I, where :, a = ;<, ;>, 8?, =. and j ≠ i.

M AN U

692

Transport term

anaerobic biomass

ACCEPTED MANUSCRIPT 42 695

3. Readily biodegradable soluble substrate (Ÿ˜s )

m87, 8, 8, 8, , , + $ + , − -j f , +$ + ,, = −f , , $ , $ , , , , mb + 8 + 8 + 8 + 8 ™šššššššššššš›ššššššššššššœ ™ššššššššššššš›šššššššššššššœ , , , , Aerobic growth of heterotrophic biomass

Anoxic growth of heterotrophic biomass

8,

$ +$ +$ + ,, −f, , , , , + 8 + 8 + , , , ™ššššššššššššššššššš›šššššššššššššššššššœ , ,

Growth of anaerobic biomas

+

*

RI PT

')

,

, ^,W + ( , + - , ) + , )87,,W − 87, * + #", $ , , , ™ššš›šššœ 6 ') , + , * , + ™ššššššššššššššš›šššššššššššššššœ ,

,

,

,

Flow in

M AN U

Transport term

4. Slowly biodegradable particulate substrate (•˜Ÿ )

, m , ^ ^,W H , , ( ) ( ) = ™š 1š −šš›š Z[ Yšššœ 1 − Z[ Y + , ) ,,W − , * + 9_, ) 7 − 7, * , + ™šššš›ššššœ , ™š š š š›š šššœ mb 6 6 Decay of Decay of ™šššššššššš›ššššššššššœ Transfer flow heterotrophic biomass

anaerobic biomass

TE D

696

,

¡¢£¤¥¡¦§¦ ¤¨ ©ª«£¬--©¢ -¬£«§®¯¥¬«© ¤£°¬ª§®¦

^ H + 9_, )87 − 87, * ™ššš›šššœ 6 H→,

,

SC

,

H→,

Flow in

, ') , , + , , * + ( , , + - , , ) − #", $ , , , ' ) * + + , ™ššššššššššššššš›šššššššššššššššœ , ,

EP

5. Soluble oxygen (Ÿ˜u )

, 8, 8, m8 , ^,W , , $ + $ + * )8 ,,W − 8 , * = −f , + X)8 − 8 + , ] šššš›š , , ™š ššššœ ™š , , , ššš›š šššœ mb 6 ™šššššššššššš›ššššššššššššœ , +8 , +8

AC C

697

¡¢£¤¥¡¦§¦ ¤¨ ©ª«£¬--©¢ -¬£«§®¯¥¬«© ¤£°¬ª§®¦

Consumption by aerobic growth of heterotrophic biomass

^ H + 9_, )8 − 8 , * ™ššš›šššœ 6 H→,

Transfer flow

Transfer term

Flow in

ACCEPTED MANUSCRIPT 43

699

6. Inert particulate (•˜± ) m , = mb

, Z Y šœ [š›š ™š ,

Decay of heterotrophic biomass

, ^,W ^ H , ) ,,W − , * + 9_, ) − , * + Z Y + [ š›ššœ ™š , , ™š 6 ššš›ššššœ ™š 6 ššš›ššššœ H→,

Decay of anaerobic biomass

Flow in

Transfer flow

7. Soluble inert organic (Ÿ˜± )

RI PT

698

, ^ m8 , ^,W H = , )8 ,,W − 8 , * + 9_, )8 − 8 , * ™ššš›šššœ ™ššš›šššœ mb 6 6 H→,

Flow in

, mg , ^,W = ± ³ mb ; Flow in

^9_ µ ; H→,

Transfer flow

701 702

SC

8. Height (depth) of the layer (²˜ )

M AN U

700

Transfer flow in

9. Recalculation of initial conditions (instant settling of solids)

Figure B1 explains the simplified approach for describing the settling of solids produced

704

by microbial growth. This calculation is used after each day of integration. Solids redistribution

705

is not considered in the ice layer.

AC C

EP

TE D

703

ACCEPTED MANUSCRIPT 44 706 707 708

RI PT

709 710 711 712

SC

713 714

M AN U

715 716 717 718

722 723 724 725 726 727 728

EP

721

AC C

720

TE D

719

Figure B1. Flow chart describing the particulates and solubles redistribution between layers for

the simulation of instant solids settlement (where ˜ = •¶, •ž, Ÿ·). Note: Excess particulates are transferred to the sludge layer without affecting the anaerobic liquid layer.

ACCEPTED MANUSCRIPT 45 729

10. Flow rates distribution

730

Input wastewater is distributed between the aerobic, ice, anaerobic, and sludge layers by using

731

the following conditions for selecting the flow rate between layers.

cd cd cd * b ∈ ¹b,W, , b,W, + Δb,W,

2. Icecap formation

cd ∈ ¹b,W, cd

b Δb

3. Icecap saturation

cd cd ¿ b ∈ ¹b,W, , beWK ‡

cd

(b) > 0

g

(b) > 0

−g

d

Input flow

^,W ∗ = 2 ¼ %L − L & ^,W ^,W∗ = h^,W ^,Wcd = 0

g

(b) > 0

^,Wd = 0 ^,W = ^,W ∗ ^,W = ^,W∗ ^,Wcd = ^,Wd∗

Transfer flow

^9_d→ ^9_ ←

Ň…

d

= −^,Wd∗ = ^,Wd∗

^9_ → cd = ^,Wd∗ − ¾2 ^9_cd← = −^,Wd∗ + ¾2 ^9_ → cd = ^,Wd∗ ^9_cd← = −^,Wd∗ ^9_ ← cd = ^,W ∗ + ¾L ^9_cd← = −^,W ∗ − ¾L

^,Wd = 0 ^,W = 0 g (b) = 0 ∗ None ^ ,W = ^,W g cd (b) < g cd cd d∗ ^,W = ^,W + ^,W ∗ cd ) and ¾L = ÂÃ · (g cd (beWK − Δb cd ) − g d ) Á9

EP

where ¾2 = Á9 ÂÃ · (g

AC C

736

cd beWK −

(b) > 0

g

TE D

4. Icecap limited

735

*

+

cd Δb,W, ,

cd cd cd ƒ b ∈ (beWK − ΔbeWK , beWK

4. Icecap melt

733 734

g

^,Wd∗ = 2 ¼ %L − L & ^,W

RI PT

cd * or b ∈ b ∈ ¹0, b,W, cd (beWK , 350ƒ

1. No icecap

732

Layer thickness

SC

Time

M AN U

Condition

ACCEPTED MANUSCRIPT

•

Mathematical model simulates biodegradation of wastewater in Arctic waste stabilization ponds. WSP model accounts for periodic operation and existence of icecap.

•

Model is calibrated and validated with data from three Arctic lagoons in Nunavut,

RI PT

•

Canada.

Model reveals significant anaerobic biodegradation of organic materials throughout

EP

TE D

M AN U

SC

winter months.

AC C

•