POLITECNICO DI MILANO

Scuola di Ingegneria Civile, Ambientale e Territoriale

Corso di Laurea Magistrale in Ingegneria per l’Ambiente e il

Territorio

IN–SITU BIOLOGICAL BIOGAS UPGRADING:

MODEL VALIDATION AND PARAMETRIC

IDENTIFICATION

Supervisor: Prof. Ing. Francesca Malpei

Co-supervisor: Ing. Viola Corbellini, PhD

Advisors: Ing. Anna Santus

Prof. Ing Gianni Ferretti

Master thesis of:

Volodymyr Maslennikov 898258

Academic Year 2019 – 2020

1

SOMMARIO

La domanda per una maggiore implementazione delle fonti rinnovabili ne sottolinea

anche un loro limite: la loro produzione è sovente discontinua. Per far fronte allo squilibrio

fra la domanda dell’energia elettrica e l’offerta fornita dalle fonti sono state sviluppate

molteplici soluzioni basate sullo stoccaggio e sulla conversione dell’energia. Uno di questi

metodi è il power-to-gas (P2G) ove si usa il surplus elettrico per operare idrolizzatori che, a

partire dall’acqua, creano idrogeno. L’idrogeno, essendo un gas dalle scarse proprietà di

stabilità e trasportabilità può essere usato per creare gas più stabili e più energeticamente

densi come metano.

La digestione anaerobica di un substrato organico come fanghi provenienti dal

trattamento delle acque municipali, produce un fango stabilizzato e biogas, una miscela

gassosa di metano, anidride carbonica ed altri gas di traccia. Sviluppando un connubio fra il

P2G e il processo della digestione anaerobica si può ottenere una procedura, chiamata bio-

upgrading che può aumentare il tenore di metano nel biogas usando l’idrogeno e l’anidride

carbonica presente nel biogas. Tale procedura viene mediata dal consorzio batterico degli

idrogenotrofi metanigeni che a partire dall’idrogeno e anidride carbonica producono

metano. L’upgrading può essere effettuato in un reattore di digestione anaerobica

preesistente oppure in un reattore speciale. Nel primo caso l’upgrading è detto in-situ, nel

secondo ex-situ

La descrizione del processo di digestione anaerobica è effettuata da svariati modelli

matematici di cui i più rilevanti sono ADM1 e AMOCO. La modellizzazione del processo di

bio-upgrading in-situ mediato da metanigeni idrogenotrofi fu effettuata da (Rosito, 2019)

mediante un modello AMOCO, modificato con equazioni dell’ADM1.

Lo scopo di questo lavoro è di validare il modello ottenuto. Ciò è stato effettuato

migliorando le equazioni del modello, valutando la risposta a diversi segnali ed effettuando

2

test di sensitività parametrica. La validità del modello come strumento di previsione è stata

testata poi con dati reali provenienti da un esperimento di upgrading in-situ. Infine, è stato

scritto un algoritmo di identificazione parametrica basato sulla trasformazione lineare

fratta.

3

ABSTRACT

The demand for a greater implementation of the renewable energy in the context of

power grid network underlines its limits: its production is often unreliable and

discontinuous. To overcome the challenge of disconnect between the demand of electric

energy and supply provided by the renewable energy sources, various solutions were

developed. One of these solutions is power-to-gas (P2G). The excess electric energy is used

to hydrolyze water to create oxygen and hydrogen. Since the hydrogen gas is characterized

by high instability and low transportability it can be used in-loco to create a more stable and

energy dense methane.

The anaerobic digestion of an organic substrate such as municipal wastewater sludge

harnesses biogas, a mixture of methane, carbon dioxide and other trace compounds. By

developing a synergy between P2G and anaerobic digestion, a procedure can be obtained

that raises the methane concentration in the biogas while removing CO2. Such procedure is

called biogas bio-upgrading and is mediated by methanogenic hydrogenotrophic microbial

communities that synthetize methane from hydrogen and carbon dioxide. The upgrading is

performed by injecting hydrogen gas either in an existing biological reactor (in -situ

procedure) or in ad-hoc vessel where no sludge is provided (Ex-situ)

The mathematical description of the anaerobic digestion process is performed by

various models. The most relevant ones are the ADM1, developed by international water

association and AMOCO, written by a research team for a European union project. The

modelling of the bio-upgrading process for the in-situ implementation was first performed

by (Rosito, 2019) by combining AMOCO model with some of ADM1’s equations.

The aim of this work is the validation of such model. The validation was performed

by testing the response to various input signals, by implementing different types of reactor

control, parametric sensitivity testing and mass balance checks. The model’s potential as a

predictive tool was then tested by comparing it against the existing data of a bio-upgrading

4

experiment. At last, a parametric identification algorithm was developed on basis of the

linear fractional transformation framework

5

TABLE OF CONTENTS

Sommario ................................................................................................................................... 1

Abstract ...................................................................................................................................... 3

List of Figures ............................................................................................................................. 8

List of Tables ............................................................................................................................ 11

Introduction ............................................................................................................................. 12

1 State of the art .................................................................................................................. 16

1.1 Anaerobic Digestion .................................................................................................. 16

1.1.1 Hydrolysis & Disintegration ............................................................................... 17

1.1.2 Acidogenesis ...................................................................................................... 17

1.1.3 Acetogenesis ...................................................................................................... 18

1.1.4 Methanogenesis................................................................................................. 18

1.1.5 Homoacetogenesis: an additional process. ....................................................... 18

1.1.6 Parameters influencing the process .................................................................. 18

1.2 Biogas Production ..................................................................................................... 21

1.2.1 Biogas composition ............................................................................................ 22

1.2.2 The Effects of the Pollutants .............................................................................. 23

1.2.3 Threshold Values for achievement of biomethane ........................................... 24

1.3 Upgrading Technologies ............................................................................................ 25

1.3.1 Physio - Chemical Upgrading ............................................................................. 25

1.3.2 Scrubbing ........................................................................................................... 25

6

1.4 Biological Upgrading.................................................................................................. 28

1.4.1 In-Situ reactor configuration .............................................................................. 29

1.4.2 Ex-situ process ................................................................................................... 30

1.5 Anaerobic Digestion Modeling .................................................................................. 30

1.5.1 ADM1 ................................................................................................................. 30

1.5.2 AMOCO .............................................................................................................. 36

2 Materials & Methods ........................................................................................................ 39

2.1 Experimentation ........................................................................................................ 39

2.2 Model Validation ....................................................................................................... 41

2.3 AMOCO model overview ........................................................................................... 41

2.3.1 Biological process overview ............................................................................... 43

2.3.2 Model control .................................................................................................... 46

2.3.3 Equations ........................................................................................................... 47

2.3.4 Parameters ......................................................................................................... 54

2.4 LFT Model .................................................................................................................. 63

2.5 AMOCO model linearization ..................................................................................... 66

3 Results and Discussion ...................................................................................................... 68

3.1 Sludge Modulation .................................................................................................... 69

3.1.1 Total Solid Output .............................................................................................. 70

3.1.2 Soluble COD output ........................................................................................... 70

3.1.3 Biomass .............................................................................................................. 71

3.1.4 Inorganic Carbon ................................................................................................ 71

3.1.5 pH ....................................................................................................................... 72

3.1.6 Gas Phase Response .......................................................................................... 72

3.1.7 Mass Balances .................................................................................................... 73

3.2 Parameter Sensitivity: VFA accumulation ................................................................. 75

7

3.3 Hydrogen injection Implementation ......................................................................... 79

3.3.1 On exceeding the stoichiometric ratio .............................................................. 80

3.4 pH Control Implementation ...................................................................................... 83

3.5 Data Fit ...................................................................................................................... 84

3.5.1 Total Solids ......................................................................................................... 84

3.5.2 Soluble COD ....................................................................................................... 85

3.5.3 Alkalinity ............................................................................................................. 85

3.5.4 Gas Phase ........................................................................................................... 86

3.6 LFT parameter Identification ..................................................................................... 88

4 Conclusions ....................................................................................................................... 90

Bibliography ............................................................................................................................. 92

Appendix A ............................................................................................................................... 95

A.1 Variables ......................................................................................................................... 95

A.2 Petersen Matrix.............................................................................................................. 96

A.3 Model’s Equations ......................................................................................................... 97

Appendix B ............................................................................................................................. 100

OpenModelica Code ........................................................................................................... 100

Note on OpenModelica version ..................................................................................... 100

Appendix C ............................................................................................................................. 134

C.1 LFTfun Function ........................................................................................................ 134

C.2 Main Script ................................................................................................................ 141

C.3 Load_Parameters script ............................................................................................ 145

C.4 EquationsToMatrix script ......................................................................................... 147

8

LIST OF FIGURES

Figure 1: CO2 concentration trend in Holocene era (IPCC, 2014) ........................................... 12

Figure 2: World energy consumption (British Petroleum Report, 2018) ................................ 13

Figure 3: summer and winter mean consumption and photovoltaic generation (PV): (a)

Lisbon summer, (b) Lisbon winter, (c) Helsinki summer, (d) Helsinki winter. (Paatero, et al.,

2007) ........................................................................................................................................ 14

Figure 4: P2G within biogas upgrading framework (Angelidaki et al, 2018) ........................... 15

Figure 5: Growth Rate as a Function of Temperature ............................................................. 19

Figure 6: pH dependence of the Methanogenic Biomass ....................................................... 20

Figure 7:Typical Values of Raw Biogas ..................................................................................... 23

Figure 8: Minimum Characteristics of the biomethane for the Italian natural gas grid (DM

2/3/18) ..................................................................................................................................... 25

Figure 9: Column Absorption Process Overview ..................................................................... 26

Figure 10: PSA Operational Scheme ........................................................................................ 27

Figure 11: Polyamide membrane gas retention scheme ......................................................... 27

Figure 12: Biogas Bio-Upgrading Scheme ................................................................................ 29

Figure 13: Process Flow Diagram - Biological Process ............................................................. 44

Figure 14: PI pH Control – the control starts at day 50. .......................................................... 47

Figure 15:Hydrogen injection rate as obtained by a PI controller ........................................... 47

Figure 16: Inorganic Carbon Speciation ................................................................................... 50

Figure 17: Un-ionized ammonia percentage as a function of pH and temperature (Speece,

1996) ........................................................................................................................................ 51

Figure 18: IpH function illustration - comparison with a pure gaussian function..................... 54

9

Figure 19: particulate COD input ............................................................................................. 55

Figure 20: Solid Matter Output ................................................................................................ 70

Figure 21: Soluble COD output – comparison with the initial conditions ............................... 70

Figure 22: Soluble COD output - up close ................................................................................ 71

Figure 23: Biomass Output Concentration .............................................................................. 71

Figure 24: Inorganic Carbon ..................................................................................................... 72

Figure 25: pH ............................................................................................................................ 72

Figure 26: gas flowrate signal in response to sludge modulation ........................................... 73

Figure 27: Gas fractions in the outflow as response to sludge modulation signal .................. 73

Figure 28: COD mass balance ................................................................................................... 74

Figure 29: Nitrogen Mass Balance ........................................................................................... 74

Figure 30: Carbon Mass Balance .............................................................................................. 75

Figure 31: Sankey Diagram for the COD flow .......................................................................... 75

Figure 32: Comparison between the real VFA output(blue) and the model output(blue)

[mgCOD/l] ................................................................................................................................ 76

Figure 33: pH signal with no inhibition .................................................................................... 77

Figure 34: VFA output with the inhibition dynamic present. .................................................. 77

Figure 35: pH with inhibition present ...................................................................................... 78

Figure 36: Inorganic Carbon stripping ..................................................................................... 78

Figure 38: Hydrogen Inflow ..................................................................................................... 79

Figure 39: Inorganic Carbon ..................................................................................................... 80

Figure 40: Biomass Washout ................................................................................................... 81

Figure 41: Output Gas Fractioning - no bypass ........................................................................ 81

Figure 42: Influent Gas Fractioning with bypass parameter ................................................... 82

Figure 43: pH signal comparison .............................................................................................. 83

Figure 44: Total Solids in the outflow - data vs model ............................................................ 84

Figure 45: Total Solids in the inflow – data vs model .............................................................. 84

Figure 46: Total Soluble COD ................................................................................................... 85

Figure 47: Alkalinity - Comparison between the model and data ........................................... 85

Figure 48: Total Gas Outflow - model vs data ......................................................................... 86

Figure 49: methane and carbon dioxide gas outflow .............................................................. 86

Figure 50: Methane and Carbon Dioxide fractioning in the output gas .................................. 87

10

Figure 51: Hydrogen Gas Outflow Comparison ....................................................................... 87

Figure 52: Error Dynamic as a function of simulation time ..................................................... 89

11

LIST OF TABLES

Table 1 Anaerobic Degradation Process (Riffat, 2013) ............................................................ 16

Table 2: Methanogen strains (DegremontTM water handbook) .............................................. 22

Table 3: Biogas impurities and their consequences (Ryckebosch et al., 2011). ..................... 24

Table 4: ADM1 components of the liquid phase ..................................................................... 32

Table 5:: ADM! Biochemical Processes .................................................................................... 33

Table 6: AMOCO variables ....................................................................................................... 38

Table 7: Model's parameters ................................................................................................... 38

Table 8: H2/CO2 ratios for each phase and the resulting hydrogen gas injection rate .......... 40

Table 9: Experimentation Data ................................................................................................ 41

Table 10: State Variables of the modified AMOCO model ...................................................... 43

Table 11: Petersen matrix ........................................................................................................ 46

Table 12: Bresso sludge characterization (Rosito, 2019) -....................................................... 55

Table 13: Input parameters ..................................................................................................... 58

Table 14: Biochemical parameters .......................................................................................... 59

Table 15: Physiochemical parameters ..................................................................................... 60

Table 16: Reactor parameters ................................................................................................. 61

Table 17: Fractioning parameters ............................................................................................ 62

Table 18: Inhibition parameters .............................................................................................. 62

Table 19: Sludge Input Signal ................................................................................................... 68

Table 20: Concentrations of the Input Species ........................................................................ 69

Table 21: Parameter Tweaking Results .................................................................................... 76

Table 22: Hydrogen bubbling Values ....................................................................................... 79

Table 23: Gas Fractioning -: Comparison between the data and the model ........................... 82

Table 24: Error as a function of time ....................................................................................... 89

INTRODUCTION

Climate change is one of the world’s most pressing challenges. Human emissions of

greenhouse gases – carbon dioxide (CO2), nitrous oxide, methane, and others – have

increased global temperatures by around 1℃ since pre-industrial times. A changing climate

has a range of potential ecological, physical and health impacts, including extreme weather

events (such as floods, droughts, storms, and heatwaves); sea-level rise; altered crop

growth; and disrupted water systems (IPCC Fifth Assessment Report, 2014).

The emission of CO2 is especially problematic since almost every human activity

produces consistent amounts of this gas. Its emission is estimated to be around 1000 kg/s

(Adnan & Ong, 2018). The entity of such emission has led to its atmospheric concentration

never seen in recent planetary history-

Figure 1: CO2 concentration trend in Holocene era (IPCC, 2014)

The main source of carbon dioxide emissions is the energy production from fossil

fuels such as oil, coal and natural gas. While the current public and academic perspective is

oriented towards finding more sustainable energy sources, it can’t be denied that the trend

Atmospheric CO2 concentration trend throughout

Introduction

13

of fossil consumption will be kept positive for the foreseeable future (BP report, 2018),

especially since more third world countries will participate in the global economy.

Figure 2: World energy consumption (British Petroleum Report, 2018)

Technologies such as Carbon Capture and Storage (CCS) may help reduce the

greenhouse gas emissions generated by fossil fuels, and nuclear energy can be a zero-

carbon alternative for electricity generation. But other, more sustainable and less risky

solutions exist: energy efficiency and renewable energy. Renewable energy is derived from

resources that are replenished naturally on a human timescale. Such resources include

biomass, geothermal heat, sunlight, water, and wind. All these sources have their strengths

and weaknesses. Some are more suited to certain locations than others. For instance some

only produce electricity intermittently, others, such as biomass, hydropower, and

geothermal, can be used as baseload generation, producing a constant, predictable supply

of electricity.

Since a prolonged storage of electric energy is costly and inefficient, a wide

application of intermittent renewable energy sources is difficult (Pataro et al, 2007).

Discrepancies between the supply provided by photovoltaic panels and the demand in two

major cities can be observed Figure C.

World Fuel Consumption - Million tonnes oil

Introduction

14

Figure 3: summer and winter mean consumption and photovoltaic generation (PV): (a) Lisbon summer, (b) Lisbon winter, (c) Helsinki summer, (d) Helsinki winter. (Paatero, et al., 2007)

An innovative approach to the matter is called Power-to-Gas (P2G). P2G uses excess

electric energy to operate an electrolyzing unit which splits water into oxygen and hydrogen

with the use of high voltage (Bunger et al., 2014). The hydrogen gas can be utilized

efficiently in power cells and in combustion units alike, but it has several drawbacks.

Although it has an excellent energetic density in regard to unit of mass, to achieve a

reasonable volumetric energetic density it has to be brought to high pressures. More than

that, since the hydrogen gas is one of the smallest molecules that exist, its prolonged

storage always means mass loss with possible ignition incidents.

These problematics prompt the conversion of hydrogen to a more stable yet energy

dense gas such as methane, so a process of electric energy conversion mediated by

hydrogen is named power to methane (P2M).

Power-to-gas systems may be deployed as adjuncts to wind parks or solar-electric

generation but are especially interesting in the context of anaerobic digestion. The

anaerobic digestion process yields a renewable resource, biogas which is a mixture of

methane and carbon dioxide. Since this mixture is unsuitable to be used as-it-is, it must be

cleaned, and the carbon dioxide concentration lowered. This biogas cleaning is called

upgrading. The coupling of anaerobic digestion with a P2G system can provide the

Introduction

15

upgrading of biogas by converting the endogenous carbon dioxide and the hydrolysis-

supplied hydrogen to obtain a biogas of higher quality.

Figure 4: P2G within biogas upgrading framework (Angelidaki et al, 2018)

Chapter 1

1 STATE OF THE ART

1.1 ANAEROBIC DIGESTION

The Anaerobic digestion is an organic matter degradation process that is carried on by

a bacterial consortium in the absence of oxidizing oxygen-based compounds such as O2 and

NO3-. It is a complex process, involving many different steps that ultimately lead to the

production of biogas:

Table 1 Anaerobic Degradation Process (Riffat, 2013)

State of the art

17

Because the biochemical redox reactions that take place in the reactor must use oxidants

that grant low net positive energy balance, or low net Gibbs energy (Speece, 1997)., that

means the microbial communities are less capable to produce large biomass sludge volume,

requiring thusly less micronutrients that are only removed due to the needs related to the

biomass synthesis (Malpei & Gardoni, 2008). The anaerobic process is mainly an

endothermic process, meaning that the heat must be supplied externally, which is the main

drawback of anaerobic digestion. As it was said before, the AD is a multistep process that

can be, however be divided in a set of processes.

1.1.1 Hydrolysis & Disintegration

The organic matter input is finely minced and the complex organic molecules such as

proteins, carbohydrates and lipids that form the substrate for the anaerobic digestion are

degraded by enzymes that are produced by hydrolyzing species of bacterial consortium.

These substances are then degraded to simpler ones such as amino acids, long chain fatty

acids (LCFA) and sugar monomers. This process is considered the limiting one due to its rate

being dependent on the substrate type: complex substrates such as lignocellulosic matter

are a strong limiter to the reaction rate. The hydrolysis is also inhibited by high

concentration of amino acids and sugars.

1.1.2 Acidogenesis

The organic monomers are further degraded to form volatile fatty acids, acetic acids

and, to some extent, CO2 and H2. These transformations are performed, by intracellular

enzymes of the acidogenic bacteria. The volatile fatty acids are mainly propionic, butyric and

valeric acids. Lactic acid and ethanol are trace byproducts that are quickly converted to

organic acids. Acetic acids and CO2, H2 are also produced. Since there is a consistent acid

production, this step has the potential to increase significatively the reactor pH, event that

can lead, in extreme cases to the process collapse. Acetogens themselves are quite

insensible to the pH, tolerating oscillations in the range between 4 – 8.5 but not so other

consortium dwellers, especially methanogens.

State of the art

18

1.1.3 Acetogenesis

This process encompasses a series of bio-reactions that lead to the degradation of

volatile fatty acids into acetic acid. The parameters of this reactions also regulate the

possible accumulations of acidity which has detrimental effects on the overall process. Since

the acetogens are themselves inhibited by high hydrogen concentrations, they require the

syntrophic action of a H2-utilising species: the hydrogen, if not removed can lead to

acetogen inhibition, which leads to VFA not being processed to acetic acid, thus leading to a

pH drop.

1.1.4 Methanogenesis

The methanogenesis is the last part of the chain of reactions that lead to the

production of biogas. Methanogenesis is carried out by a group of Archaea organisms

collectively known as methanogens. There are two groups of methanogens. One group is

termed acetoclastic (or acetotrophic) methanogens, who split acetate into methane and

CO2. The second group is termed hydrogen-utilizing (or hydrogenotrophic) methanogens,

who use hydrogen as the electron donor and CO2 as the electron acceptor to produce

methane. These hydrogenotrophic bacteria utilize the H2 which the acetogens produce. H2

uptake by the methanogenic community is very efficient, having an affinity of parts per

million, which ensures very low hydrogen partial pressure (Murphy & Thamisoroj, 2003).

The Methane, being non-polar, has thus low solubility inside the water, passing almost

integrally into the gas phase.

1.1.5 Homoacetogenesis: an additional process.

The homoacetogenesis is a process that produces acetic acid and is mediated by a

community of hydrogen-consuming acetogens. Their operativity helps to keep the hydrogen

levels very low, although their H2 uptake is outcompeted by the methanogens.

1.1.6 Parameters influencing the process

The process of digestion is influenced by a multitude of parameters such as type of

reactor, Temperature, pH, presence inhibiting substances.

State of the art

19

• Reactor Typology: the most common types of reactors are the CSTR, PF and batch

reactors. The CSTR is the most common in the wastewater sludge processing, while

batch reactors see an extensive use in the putrescible organic waste stabilization.

The plug flow reactors are advantageous while dealing compounds that are toxic for

the biomass fauna (Speece, 1997) due to the fact that the PF reactors, while not

having the same mixing efficiency of the CSTR ones, are nonetheless characterized

by lower HRT at the same efficiency level, meaning that the contact time between

bacteria and toxicant.

• Temperature: the temperature influences the bacterial growth rate. From

temperature adaptation perspective there are three main groups of bacteria:

Psychrophiles, that operate at low temperatures, mesophiles that operate at

intermediate temperatures and thermophiles, that operate at high temperatures

Their respective growth rates as a function of temperature are shown in Figure 5.

Figure 5: Growth Rate as a Function of Temperature

• pH: The pH influences the metabolic processes of the bacterial population as well as

the shift of equilibrium of the organic acids and ammonia alike.

The pH of the process is dictated by the most sensitive consortium

Temperature [C°]

State of the art

20

group: the methanogens. Their activity is the highest when the pH is within the 6.8 –

8.2 range (Speece, 1996).

Figure 6: pH dependence of the Methanogenic Biomass

• Ammonia. The ammonium cation is an important nitrogen source for the biomass,

but the ammonia equilibrium, in particular conditions (namely, low pH) can favor the

unionized phase that leads to metabolic inhibition. Its control is important to avoid

the biomass intoxication.

• Inorganic carbon availability: The inorganic carbon is a multi-faceted aspect of the

wastewater treatment and anaerobic digestion alike. While on one side it serves to

increase the alkalinity, on the other side it acts as a reductant for the biochemical

reactions and carbon source for the hydrogenotrophic methanogens.

• Alkalinity: The alkalinity of a water is a measure of its capacity to

neutralize acids and is usually measured as concentration of weak acids, either in

terms of protons liberated by their dissociation or in CaCO3 equivalence. In the

anaerobic digestion there are two chemical species that cause the lion’s share of the

pH shift:

State of the art

21

o the CO2 that is the byproduct of bacterial respiration. It Increases the acidity

due to the release of a proton from an intermediate, carbonic acid of the

following chain reaction, that leads to the hydrogencarbonate production

from CO2:

{𝐶𝑂2(𝑎𝑞) + 𝐻2𝑂(𝑙) ↔ 𝐻2𝐶𝑂3(𝑎𝑞)

𝐻2𝐶𝑂3(𝑎𝑞) + 𝐻2𝑂(𝑙) ↔ 𝐻𝐶𝑂3− + 𝐻3𝑂(𝑎𝑞)

+

o VFA: anions and their salts can act as buffer against added acidity but not

from VFA itself (Speece 1996). This means that their accumulation can lead

to the pH drop

• Organic Loading Rate: it’s a measure of substrate flow that insists on reactor’s

biomass. Since the fermenting bacteria are sensitive to the volume and type of

feedstocks injected into a digester, its value is important to determine situations of

starvation or substrate-induced inhibition.

• HRT/SRT: the hydraulic residence time is the mean time of dissolved compound’s

permanence in the reactor, as their destiny is that of the water flow. A high enough

value in a situation of a single stage CSTR reactor is paramount for the feasibility of

the process. the solids retention (SRT) time is the mean permanence time of a solid

particle in the reactor. In the case of a CSTR reactor with no sludge recirculation or

retention it is coincident with HRT.

1.2 BIOGAS PRODUCTION

The main digestible components of wastewater sludge are carbohydrates, proteins

and fats and their concentration in the sludge determines the entity of methanation.

Species of only two genera, Methanosarcina and Methanothrix can produce methane from

State of the art

22

acetic acid and approximately 70% of methane produced comes from acetoclastic pathway,

with the remainder being the hydrogenotrophic way (Murphy, T. Thamisisroj, 2003).

Species Substrate

Methanobacterium H2/CO2

Methanobrevibacter H2/CO2

Methanococcus H2/CO2

Methanosarcina H2/CO2/acetate

Methanothrix acetate

Table 2: Methanogen strains (DegremontTM water handbook)

The quantity of methane that can be predicted from the well-known stoichiometric formula

(Buswell and Hatfield, 1936), where the parameters are dependent upon the sludge

characterization.

𝐶𝑛𝐻𝑎𝑂𝑏 + (𝑛 −𝑎

4+𝑏

2)𝐻2𝑂 → (

𝑛

2+𝑎

8−𝑏

4)𝐶𝐻4 + (

𝑛

2−𝑎

8+𝑏

4)𝐶𝑂2 1.1

1.2.1 Biogas composition

The nature of the raw materials and the operational conditions used during

anaerobic digestion, determine the chemical composition of the biogas. Biogases with

varying degrees of purity and composition can be produced at different facilities such as

sewage treatment plants (sludge fermentation stage), landfills, sites with industrial

processing industry and at digestion plants for agricultural organic waste, both mesophilic

(35 °C) and thermophilic (55 °C) (Ryckebosch et al., 2011).

Raw biogas consists mainly of methane (50–80%) and carbon dioxide (30–50%).

Trace amounts of other components such as water (H2O, 5–10%), hydrogen sulfide (H2S,

0.005–2%), siloxanes (0–0.02%), halogenated hydrocarbons (VOC, < 0.6%), ammonia (NH3,

<1%), oxygen (O2, 0–1%), carbon monoxide (CO, <0.6%) and nitrogen (N2, 0–2%) can be

present. The values for the typical raw biogas composition are reported in the table below

(Chen et al., 2015)

State of the art

23

Figure 7:Typical Values of Raw Biogas

1.2.2 The Effects of the Pollutants

The composition of the biogas strongly impacts the feasibility of its utilization: The

high CO2 content, lowers the gas calorific value while the water reacts with hydrogen sulfide

and ammonia to form acids (Angelidaki et al., 2013). More than that, some of the trace

pollutants create a separate set of problems: the siloxanes are able to create occluding

residues, while hydrocarbons can create corrosion problems to the engines when ignited

(Ryckebosch et al., 2011). The common issues with the biogas pollutants are reported in

Table 3.

Impurity Possible Impact

Water

Corrosion in compressors, gas storage tanks and engines due to reaction with

H2S, NH3 and CO2 to form acids

Accumulation of water in pipes

Condensation and/or freezing due to high pressure

Dust Clogging due to deposition in compressors, gas storage tanks

𝐻2S Corrosion in compressors, gas storage tanks and engines

Toxic concentrations of H2S (> 5 𝑐𝑚3 𝑚3⁄ ) remain in the biogas

SO2 and SO3 are formed due to combustion, which are more toxic than H2S and

cause corrosion with water

𝐶𝑂2 Low calorific value

Siloxanes Formation of SiO2 and microcrystalline quartz due to combustion; deposition

at spark plugs, valves and cylinder heads abrading the surface

State of the art

24

Hydrocarbons Corrosion in engines due to combustion

𝑁𝐻3 Corrosion when dissolved in water

O2/air Explosive mixtures due to high concentrations of O2 in biogas

𝐶𝑙− Corrosion in combustion engines

𝐹− Corrosion in combustion engines

Table 3: Biogas impurities and their consequences (Ryckebosch et al., 2011).

1.2.3 Threshold Values for achievement of biomethane

Because of these issues the Biogas has to undergo several processes that remove the

problematic compounds. The consistent removal of these compounds, especially with

regards to CO2 produces a fuel-grade substance, the Biomethane. The Italian normative (DM

2/3/18) specifies the threshold values for a biogas to be considered biomethane, as

reported in the table below.

Parameter Acceptable values Unit of measure

Higher heating

value

34.94 ÷ 45.28 MJ/Sm3

Wobbe Index 47.31 ÷ 52.33 MJ/Sm3

Relative Density 0.5548 ÷ 0.8 -

Dew Point of water

(at 7000 kPa)

≤ -5 °C

Dew Point of

Hydrocarbons

≤ 0 °C

Tmax <50 °C

Tmin >3 °C

Oxygen ≤0.6 % mol

Carbon dioxide ≤3 % mol

Hydrogen Sulfide ≤6.6 mg/Sm3

Hydrogen ≤0.5 % mol

Carbon monoxide ≤0.1 %mol

Sulfur from

mercaptans

≤15.5 mg/Sm3

Total sulfur ≤150 mg/Sm3

Mercury ≤1 μg/Sm3

State of the art

25

Chlorine <1 mg/Sm3

Fluorine <3 mg/Sm3

Silicon ≤5 ppm

Figure 8: Minimum Characteristics of the biomethane for the Italian natural gas grid (DM 2/3/18)

1.3 UPGRADING TECHNOLOGIES

The biogas is treatable a number of different technologies. They can be categorized

into two main groups: the physio-chemical and biological technologies. The physio -

chemical methods involve consolidated technologies such as membrane filtering and

column scrubbing. They allow both the cleaning and the upgrading to be obtained but they

are costly and require installation of new operational units. The Biological technologies on

the other hand can be installed in pre-existing facilities and allow the utilization of biogas’

CO2 for the methane synthesis but do not allow the gas cleanup from other compounds.

1.3.1 Physio - Chemical Upgrading

The physio-chemical methods aim to remove the carbon dioxide by the use of five

mainstream technologies, involving processes of adsorption, absorption and membrane

separation. These technologies are capable of achieving methane recovery as high as 96%.

1.3.2 Scrubbing

Scrubbing is one of the most common processes in the industrial chemistry and

sanitation engineering alike and it exploits the fact that CO2 and other trace pollutants such

as H2S have higher water solubility than The mass transfer occurs inside a column reactor

where the gas phase and the solvent, water are being flown counter-current. The water is

sprayed out of top positioned nozzles, usually on top of the packing material, while gas is

being blown from the bottom. The packing material acts as surface – increaser between the

two phases. The water is a close circuit vector that has to be opportunely treated to remove

the absorbed pollutants.

Scrubbing by using organic solvents

By using solvents other than water the difference of solubility between the pollutant

and the methane is emphasized. Organic solvents such as methanol, and dimethyl ethers

are used. Proprietary solvents such as SelexolTM have solubility three times higher than

water

State of the art

26

Chemical Scrubbing

The chemical scrubbing involves the use of a solvent that reacts with the pollutant of

interest, creating chemical bonds. In the case the CO2, aqueous amine solutions are

commonly used. These solutions have the advantage that the H2S is also absorbed. The

solvent has to be regenerated in a stripping process where it is subjected to temperatures

as high as 160 °C and pressures in the range of 1-2 bar to disrupt the bonds between the

CO2 and amine solution.

Figure 9: Column Absorption Process Overview

Pressure swing adsorption

The adsorption is a process of formation of the weak Wan-der-Walls bonds between

the target gaseous or aqueous phase compound and a solid material, characterized by a

high surface per volume unit. The most common materials are activated carbon, zeolites

and carbon molecular sieve (Augelletti et al., 2017) Pressurized gasses tend to be attracted

to solid surfaces so the biogas flow is kept at high pressures during the adsorption phase (4-

10 bar), and at low pressure during the subsequent desorption phase. The adsorption

material selectively retains N2 O2 CO2 H2O H2S while methane is able to pass through

unscathed. When the adsorbent column is exhausted the adsorption, cycle ends. The main

gas flow is diverted from the said column and a clean, hot, low pressure gas is flown,

desorbing thusly the pollutants.

State of the art

27

Figure 10: PSA Operational Scheme

Membrane separation

The separation of gaseous mixture can be performed by a solid membrane that acts

as a diffusion rate classifier, allowing faster rates of membrane crossing for certain

molecules and slower for others. The process occurs between two flows that are separated

by a membrane: the feed flow, which has a high pressure and the permeate flow, at a lower

pressure. The pressure gradient thusly insured is the force driving the whole process. The

total membrane transfer rate depends upon the diffusion rate and sorption rate. The

diffusion is dependent upon the structure of the molecule, its size, its interaction with the

membrane medium, while the sorption rate, that is, the passage from the gaseous to the

solid phase is a function condensability and, to a lesser extent the rate of passage from the

bulk zone to the Prandtl limit stratum (Bauer et al., 2013).

In the context of biogas upgrading the methane has the slowest overall transfer rate with

the membranes made of cellulose acetate and polyamide (Baker et al, 2012), while CO2 the

fastest. It follows thus that the retentate is mainly CH4, while the permeate is mainly CO2.

Figure 11: Polyamide membrane gas retention scheme

State of the art

28

Cryogenic Separation

The cryogenic process uses the different molecular condensation rate as a driving force of

the separation process. The Gas mixture temperature is lowered to as low as -110°C

allowing the condensation of methane and leaving unscathed other mixture species.

1.4 BIOLOGICAL UPGRADING

The use of microbial communities to operate the biogas upgrading has several

advantages in respect to the physio-chemical upgrading. The CO2 is converted rather than

removed, into other energy form, methane and the whole process is performed at

atmospheric-like conditions, which lowers operational costs. More than that, the use of

physio – chemical processes involve additional costs due to external equipment such as

pumps, compressors and membranes.

This method involves the use of hydrogenotrophic methanogens that catalyze the

Sabatier reaction (eq. 1.2)

4𝐻2(𝑎𝑞) + 𝐶𝑂2(𝑎𝑞) → 𝐶𝐻4(𝑎𝑞) + 𝐻2𝑂(𝑙) 1.2

The hydrogen gas is supplied exogenously and must be produced in a way that makes the

whole process economically sound, e.g. it cannot be produced by using fossil fuels, which

would comport a mere transfer of energy in best case. Due to the fact that hydrogen gas is

one of the smallest molecules, it has a high diffusion rate as well as very high volumetric

density. That means that it cannot be stored efficiently in the containers made of common

materials, and prolonged storage always implies a loss of matter. Thai means it has to be

produced in loco. The most common way is the hydrolyzation of water to produce hydrogen

gas and oxygen:

𝐻2𝑂(𝑣) → 𝐻2(𝑔) + 𝑂2(𝑔) 1.3

In this way, excess electric energy can be used, either deriving from renewable sources, such

as wind and solar (photovoltaic and photothermic) or as an electric grid surplus energy

available at the times of low demand (think of a nuclear plant, that has a very low power

output variation capability throughout the day) This methodology is called power to gas

State of the art

29

(P2G) and is an attractive perspective especially when integrated with solar and wind

technologies (Kougias et al., 2017)

Figure 12: Biogas Bio-Upgrading Scheme

In regard to the reactor implementation, the chemoautotrophic methanation is

feasible in three ways: the in-situ, the in-situ and hybrid configurations.

1.4.1 In-Situ reactor configuration

In this configuration the hydrogen gas injection is implemented in an existing sludge

digester. The CO2 is obtained from the organic matter that is fed in the reactor. This

configuration is advantageous because it can be implemented with minor modifications in a

vast array of existing digesters, both high capacity industrial plants and local factories (Li, et

al., 2020). If the parameters are carefully monitored the upgrading process can produce a

biogas as 99% pure (Angelidaki et al. 2018). Significant problems, however, were underlined

by (Bassani et al., 2016; Luo & Angelidaki, 2012; Wang et al., 2013). The main concerns

associated with the process are VFA when the H2 concentration is high and a sustained

consumption of CO2, which can lead to pH increase.

State of the art

30

1.4.2 Ex-situ process

The ex-situ upgrading is carried in an ad-hoc vessel where a sludge containing only

hydrogenotrophic methanogens, that require only gaseous substratum and micronutrients,

is present. The carbon dioxide is supplied by the injection of the to-be-upgraded biogas, and

hydrogen is supplied exogenously. This process is more easily controllable, has less

dependence on substrate and there is no risk of VFA accumulation, but it requires higher

investment costs, as a reactor has to be built from scratch. The main issue regarding the

process is the low gas-liquid transfer phase which limits the rate of substratum transfer to

the liquid phase, so an adequate mixing speed has to be identified.

1.5 ANAEROBIC DIGESTION MODELING

The mathematical models can underline important aspects of the anaerobic digestion

process, such as inhibition pathways and control optimization (Lovato t al., 2017). More

than that, mathematical models minimize the experimental effort, risk and cost (Angelidaki

et al., 1999). In the field of anaerobic digestion there are two main models. Although these

models do not encompass the hydrogenotrophic archaea – driven methanogenesis, some

important steps were undertaken by Lovato et al. 2017 and Rosito, 2019 to implement this

dynamic.

1.5.1 ADM1

The IWA Anaerobic Digestion Modelling Task Group constructed a generalized

Anaerobic Digestion model (Batstone et al., 2002) called ADM1. The ultimate goal was to

support the increasing application of the anaerobic digestion process as a sustainable way

to treat wastes and produce renewable energy (Batstone et al., 2002).

The model is based upon the mass conservation principle in a CSTR reactor, which is

expressed by a set of differential equations. The biological processes of the ADM1 model are

extensive and include disintegration, hydrolysis, acidogenesis and methanogenesis,

although some minor dynamics such as formate uptake are excluded. In total, it accounts

for 19 biochemical reactions associated with 7 bacterial populations. The kinetics is

structured according to a Monod function of the substrate and takes into consideration pH,

hydrogen and ammonia inhibition terms. The model also accounts for physio-chemical

State of the art

31

reactions: liquid–gas transfer, acid–base reactions and gas phase equations. The gas mass

balance is described by three differential equations, while the ionic interactions are either

algebraic or differential.

1.5.1.1 Liquid phase equations

The sludge components can be divided into soluble and particulate components. In

total, there are 24 components in the liquid phase, 12 soluble and 12 particulate

components (Table 4: ADM1 components of the liquid phase

State of the art

32

), combining in 19 different processes (j=1–19) as shown in Table 5.

Soluble Compounds Particulate compounds

Saq,1 Ssu Monosaccharides Xaq,13 Xc Composites

Saq,2 Saa Amino acids Xaq,14 Xch Carbohydrates

Saq,3 Sfa Long-chain fatty

acids

Xaq,15 Xpr Proteins

Saq,4 Sva Total valerate Xaq,16 Xli Lipids

Saq,5 Sbu Total butyrate Xaq,17 Xsu Sugar degraders

Saq,6 Spro Total propionate Xaq,18 Xaa Amino acid degraders

Saq,7 Sac Total acetate Xaq,19 Xfa Long-chain fatty acid

degraders

Saq,8 Sh2 Hydrogen Xaq,20 Xc4 Valerate and butyrate

degraders

Saq,9 Sch4 Methane Xaq,21 Xpro Propionate degraders

Saq,10 SIC Inorganic carbon Xaq,22 Xac Acetate degraders

Saq,11 SIN Inorganic nitrogen Xaq,23 Xh2 Hydrogen degraders

Saq,12 SI Soluble inerts Xaq,24 XI Particulate inerts

Table 4: ADM1 components of the liquid phase

State of the art

33

process # Process Description

j1 Disintegration

j2 Hydrolysis of carbohydrates

j3 Hydrolysis of proteins

j4 Hydrolysis of lipids

j5 Uptake of sugars

j6 Uptake of amino acids

j7 Uptake of LCFA

j8 Uptake of valerate

j9 Uptake of butyrate

j10 Uptake of propionate

j11 Uptake of acetate

j12 Uptake of hydrogen

j13 Decay of sugar degraders

j14 Decay of amino acid degraders

j15 Decay of LCFA degraders

j16 Decay of valerate and butyrate

degraders

j17 Decay of propionate degraders

j18 Decay of acetate degraders

j19 Decay of hydrogen degraders

Table 5:: ADM! Biochemical Processes

The differential equations for the dissolved (eq. 1.4) and particulate matter (eq. 1.5) have

the following form:

dSaq,i

dt=Q

V(Sin,i − Saq,i) +∑ρjνi,j

19

j=1

i = 1 − 12 1.4

dXaq,i

dt=Q

V(Xin,i − Xaq,i) +∑ρjνi,j

19

j=1

i = 13 − 24 1.5

Where Q is the sludge inflow and the V is the volume of the CSTR reactor

State of the art

34

1.5.1.2 Acid–base equations

The acid-base equations are modelized by the implementation of carbonate,

ammonia, acetate, butyrate, valerate and propionate equilibria. If strong bases or acids are

present in the feed influent or if hydrogencarbonate is added, the model allows them to be

accounted through the cation and anion components, Scat and San.

𝑑𝑆𝑎𝑞,𝑖

𝑑𝑡=𝑄

𝑉(𝑆𝑖𝑛,𝑖 − 𝑆𝑎𝑞,𝑖); 𝑖 = 25 − 26 1.6

The ADM1’s pH is a delicate question because it is desirable to have a pure ODE

system, but since the implementation of the differential equations to express either the

hydronium or hydroxide mass balance is a difficult endeavor, the pH is usually expressed in

algebraic form. The equation 1.7 shows the ionic balance as expressed by (Batstone et al.,

2002) in the original model version, while the equation 1.8 exhibits a more sophisticated

expression as thought by (Rosen & Jeppesen, 2006). Attempts to make the ADM1 a purely

ODE system were made by (Thamsiriroj and Murphy, 2011), as expressed by the equation

1.9, which is obtained by deriving the ionic balance equation with respect to time.

Scat+ + Snh4+ + SH+ − Shco3− −Sac −

64−Spr −

112−Sbu−160

−Sva−208

− SOH− − San− = 0 1.7

{

𝜃 = 𝑆𝑐𝑎𝑡+ + 𝑆𝑛ℎ4+ − 𝑆ℎ𝑐𝑜3− −𝑆𝑎𝑐 −

64−𝑆𝑝𝑟 −

112−𝑆𝑏𝑢−160

−𝑆𝑣𝑎−208

− 𝑆𝑎𝑛−

𝑆𝐻+ =−𝜃 + √𝜃2 − 4𝑘𝑤

2

1.8

State of the art

35

dSHdt=

dSandt +

ka,IN(ka,IN + SH)

dSINdt

+ ka,IC

ka,IC + SH

dSICdt

+ka,ac

64(ka,ac + SH)

dSacdt

+ka,pro

112(ka,pro + SH)

dSprodt

+

+ka,bu

160(ka,bu + SH)

dSbudt

+ka,va

208(ka,va + SH)

dSvadt

−dSINdt

−dScatdt

1 +ka,INSIN

(ka,IN+SH)2 +

ka,ICSIC

(ka,IC+SH)2 +

ka,acSac

64(ka,ac+SH)2 +

+ka,proSpro

112(ka,pro+SH)2 +

ka,buSbu

160(ka,bu+SH)2 +

ka,vaSva

208(ka,va+SH)2 +

kw(SH)

2

1.9

1.5.1.3 Gas phase equations

Biogas considered in the ADM1 model is composed of methane, carbon dioxide and

hydrogen. The transformation of biogas from liquid to gaseous form is modelized by a Fick

mass transfer rate (equation 1.10)

𝑟𝑇,𝑖 = 𝑘𝑙𝑎,𝑖(𝑆𝑔𝑎𝑠,𝑖

𝐶𝑂𝐷𝑒𝑞,𝑖 − 𝑘𝐻,𝑖𝑝𝑖), 𝑖 = 𝐶𝑂2, 𝐶𝐻3, 𝐻2

1.10

Where CODeq quantifies the COD equivalence of a gas mole. The gas mass balance of a single

gaseous species (eq. 1.11 ) is influenced by the transfer from the aqueous phase and the

removal by the headspace gas outflow (eq. 1.12)

𝑑𝑆𝑔𝑎𝑠,𝑖

𝑑𝑡= −𝑞𝑔𝑎𝑠 ∗

𝑆𝐶𝐻𝑔𝑎𝑠,𝑖

𝑉𝑔𝑎𝑠+ 𝜌𝑇,𝑖

𝑉𝑙𝑖𝑞

𝑉𝑔𝑎𝑠

1.11

𝑞𝑔𝑎𝑠 =𝑅𝑇

(𝑝𝑔𝑎𝑠 − 𝑝𝑔𝑎𝑠,𝐻2𝑂)𝑉𝑙𝑖𝑞(𝜌𝑇,𝐻2 + 𝜌𝑇,𝐶𝐻4 + 𝜌𝑇,𝐶02)

1.12

1.5.1.4 Model Development & Applicability

The model, as presented by Batstone has been subjected by an extensive overhaul:

Blumensaat & Keller, 2005 reviewed it by correcting inconsistencies in carbon and nitrogen

mass balances while introducing equations to monitor the phosphate, TKN and alkalinity.

(Lubken et al.,2007) introduced energy balance. (Rosen et al., 2006) extended the model,

obtaining a DAE with a total of 32 differential equations and three algebraic equations.

State of the art

36

The model proved to be an excellent tool for the academic environment because of its

detailed description of the anaerobic digestion process (Thamsiriroj and Murphy, 2011). The

model inherent complexity however, limits its full-scale commercial application, either

because the full parametric set identification is not cost effective, or the modulization of

non-lab conditions is hard to obtain.

1.5.2 AMOCO

AMOCO (Bernard et al., 2001) was developed in the context of a European union

project finalized to create a simple monitoring and control tool the anaerobic digestion

process.

The model describes only two biochemical processes, acidogenesis and

methanogenesis, based upon mass balances. In the first step, acidogenesis, the acidogens

(X1) consume the organic substrate (S1) while producing CO2 (C) and VFA (S2), while in the

second step the methanogens consume the VFA to produce CO2 and CH4. The microbial

growth is expressed by the Monod dynamic in the case of the acidogens, and Haldane

dynamic in the case of the methanogens, to account for the VFA induced inhibition. The

versatility of the model allows to account for the sludge recirculation degree by the

implementation of the parameter α (0, fixed-bed; 1, ideal CSTR). The original model is thusly

formed by six ODE equations accounting for substratum, biomass carbon dioxide and

alkalinity and five algebraic equations, expressing the methane and carbon dioxide flowrate,

pH and CO2 partial pressure.

The state equations are expressed by the equations 1.13-1.18 and the ancillary algebraic

equations are expressed by equations 1.19-1.23

𝑑𝑋1𝑑𝑡

= 𝛼𝐷𝑋1 + 𝜇1𝑆1

𝑆1 + 𝐾𝑆,1 1.13

𝑑𝑋2𝑑𝑡

= 𝛼𝐷𝑋2 + 𝜇2𝑆2

𝑆2 + 𝐾𝑆,2 +𝑆2𝐾𝐼,2

1.14

𝑑𝑆1𝑑𝑡= 𝐷(𝑆1,𝑖𝑛 − 𝑆1) − 𝑘1𝜇1

𝑆1𝑆1 + 𝐾𝑆,1

1.15

State of the art

37

𝑑𝑆2𝑑𝑡= 𝐷(𝑆2,𝑖𝑛 − 𝑆2) + 𝑘2𝜇1

𝑆1𝑆1 + 𝐾𝑆,1

− 𝑘3𝜇2𝑆2

𝑆2 + 𝐾𝑆,2 +𝑆2𝐾𝐼,2

1.16

𝑑𝑍

𝑑𝑡= 𝐷(𝑍𝑖𝑛 − 𝑍)

1.17

𝑑𝐶

𝑑𝑡= 𝐷(𝐶𝑖𝑛 − 𝐶) + 𝑘4𝜇1

𝑆1𝑆1 + 𝐾𝑆,1

+ 𝑘5𝜇2𝑆2

𝑆2 + 𝐾𝑆,2 +𝑆2𝐾𝐼,2

1.18

𝑞𝑐 = 𝑘𝑙𝑎(𝐶 + 𝑆2 − 𝑍 − 𝐾𝐻𝑃𝑐) 1.19

𝑃𝐶 =−𝜑 ± √𝜑2 − 4𝐾𝐻𝑃𝑇(𝐶 + 𝑆2 − 𝑍)

2𝐾𝐻 1.20

𝜑 = 𝐶 + 𝑆2 − 𝑍 + 𝐾𝐻𝑃𝑇 +𝑘6𝑘𝑙𝑎(𝜇2

𝑆2

𝑆2 + 𝐾𝑆,2 +𝑆2𝐾𝐼,2

) 1.21

𝑞𝑚 = 𝑘6𝜇2𝑆2

𝑆2 + 𝐾𝑆,2 +𝑆2𝐾𝐼,2

1.22

𝑝𝐻 = −𝑙𝑜𝑔(𝐾𝑏𝐶 − 𝑍 + 𝑆2𝑍 − 𝑆2

) 1.23

State of the art

38

The Table 6 and Table 7 showcase the model’s variables and parameters respectively.

Variable Description Measure

X1 Acidogens kgCOD/m3

X2 Methanogens kgCOD/m3

S1 Organic Substrate kgCOD/m3

S2 VFA kgCOD/m3

Z Alkalinity M

C CO2 M

qC Gaseous CO2 flowrate m3/s

qM Gaseous CH4 flowrate m3/s

PC Carbon Dioxide Partial Pressure /

Table 6: AMOCO variables

Parameter Description Measure

D Sludge Overturn Rate (~ HRT-1) 1/d

α Sludge Recirculation Entity /

μ1 Maximum growth rate of the substrate degraders 1/d

μ2 Maximum growth rate of the VFA degraders 1/d

Ks,1 Half-saturation constant for organic substrate uptake kgCOD/m3

Ks,2 Half-saturation constant for VFA uptake kgCOD/m3

KI,2 VFA-induced inhibition constant kgCOD/m3

KH Henry constant

Table 7: Model's parameters

AMOCO has consistent simplifications in regards of ADM1. The particulate is

considered to be completely hydrolyzed and the pH is estimated roughly on the

hydrogencarbonate equilibrium. More than that, relevant ammonia, pH and H2 inhibitions

are neglected. While the model’s complexity doesn’t allow it’s use for designing and

modeling, the reductionistic approach allows however a brisk parametric identification and

implementation in a vast array of situations where process monitoring and control is

needed.

Chapter 2

2 MATERIALS & METHODS

The work of this thesis aims to expand the previous work of (Rosito, 2019), whose goal was

to perform a biogas upgrading experimentation and develop a mathematical model that

could describe such a process. The continuation of this endeavor is obtained by validating,

modifying, and expanding the model previously developed. The validation procedure is

described in Chapter 3, while in this chapter the previous experimental work is presented

and the model is presented.

2.1 EXPERIMENTATION

The aim of the experiment performed by (Rosito, 2019) was to study the in-situ

biological upgrading process applied to the anaerobic digestion of wastewater sludge. It was

performed within two lab-scale reactors at the A. Rozzi laboratory of the Department of Civil

and Environmental Engineering at Politecnico di Milano – Cremona Campus.

The CSTR reactors were fed with sludge coming from Bresso wastewater treatment plant

and the upgrading was obtained by hydrogen gas injection in the reactor in a succession of

seven increasing steps during a period of 140 days. The hydrogen inflow 𝑞,ℎ2,𝑖 for each the

upgrading phase i was defined on basis of gaseous carbon dioxide reactor outflow in the

pre-upgrading phase 𝑞𝑐𝑜2,𝑝𝑟𝑒 (phase where the reactors were let to run for 20 days without

hydrogen injection). The hydrogen gas inflow was computed as per eq. 2.1:

𝑞,ℎ2,𝑖 = (𝐻2𝐶𝑂2

)𝑖

𝑞𝑐𝑜2,𝑝𝑟𝑒 2.1

Materials & Methods

40

With (𝐻2

𝐶𝑂2)𝑖 being the ratio between the available CO2 and the supplied H2 with 4 being the

stoichiometric value for the methane synthesis. The values for each reactor are different

because the carbon dioxide outflow is different. The results of the calculation of the 𝒒,𝒉𝟐,𝒊

values for each reactor and each phase are reported in Table 8.

Starting day 0 20 41 49 58 67 94 107

Rea

cto

r 1

(𝑯𝟐𝑪𝑶𝟐

)𝒊

Pre-H2 0.76 1.01 1.9 2.7 4.36 6.01 7

𝒒,𝑯𝟐,𝒊[𝒎𝒍 𝒅⁄ ] Pre-H2 396 795 1580 1995 3438 4987 5440

Rea

cto

r 2

(𝑯𝟐𝑪𝑶𝟐

)𝒊

Pre-H2 0.8 1 1.9 2.6 4.4 6 7

𝒒,𝑯𝟐,𝒊[𝒎𝒍 𝒅⁄ ] Pre-H2 596 793 1496 2028 3438 5486 5486

Table 8: H2/CO2 ratios for each phase and the resulting hydrogen gas injection rate

A vast array of reactor and sludge parameters were obtained during the

experimentation, resulting in the creation of a time-dependent dataset that covers all the

140 days of experimentation. Its summary with regards to the identification of key sludge

components and parameters are reported in Table 9.

MEAN STDEV DESCRIPTION UNIT I/O

𝑪𝑶𝑫𝑿,𝒕𝒐𝒕 28.07 4.82 Total particulate

CDD [𝑔𝐶𝑂𝐷/𝑙] input

𝑨𝒍𝒌 1600.00 0.00 Alkalinity [𝑚𝑔𝐶𝑎𝐶𝑂3/𝑙] input

𝑻 36.89 1.33 Outflow sludge

temperature [𝐾] output

𝒑𝑯 7.34 0.13 output

𝑪𝑶𝑫𝑺,𝒕𝒐𝒕 3265.33 467.14 Total soluble COD [𝑚𝑔𝐶𝑂𝐷/𝑙] output

𝑽𝑭𝑨 812.43 640.72 [𝑚𝑔𝐶𝑂𝐷/𝑙] output

Materials & Methods

41

𝑪𝑶𝑫𝑿,𝒕𝒐𝒕 18.63 2.52 Total particulate

COD [𝑚𝑔𝐶𝑂𝐷/𝑙] output

𝑨𝒍𝒌 4002.20 297.80 Alkalinity [𝑚𝑔𝐶𝑎𝐶𝑂3/𝑙] output

𝒒𝒈𝒂𝒔 4735.14 1051.37 Gas outflow [𝑚𝑙/𝑑] output

𝑪𝑶𝟐 18.68 4.35 Carbon dioxide

fraction % output

𝑪𝑯𝟒 72.16 3.43 Methane fraction % output

𝑯𝟐 3.20 2.30 Hydrogen fraction % output

Table 9: Experimentation Data

2.2 MODEL VALIDATION

The model developed by (Rosito, 2019) is a crossbreeding between AMOCO and ADM1

models. It was implemented in OpenModelica environment. The validation, as described in

chapter 3, was performed by modulation of input, control variables tweaking and

parametric sensitivity tests. Mass balances with regards to COD, carbon and nitrogen were

carried out to test the model’s consistency. The comparison between the model and the

dataset was performed at last.

The model’s parametric identification was executed in MATLAB with the aid of an

optimization and linearization algorithm, likes of which will be discussed further on.

2.3 AMOCO MODEL OVERVIEW

The AMOCO, as described in the previous chapter, is mainly a control and monitoring tool

and as such is ill suited to be a simulating and descriptive tool. The ADM1 is on the other

hand excessively verbose and resource intensive when it comes to the dynamic simulation.

To overcome such limitations, especially with regards to the inclusion of the

hydrogenotrophic methanogen dynamic, the extension of the model was developed by

(Rosito, 2018).

The model’s foundations are the AMOCO equations and the extension is achieved by

implementing some of the ADM1 equations and parameters alike. The standard AMOCO

Materials & Methods

42

equations describe only the acidogenesis and the methanogenesis. The modified AMOCO

model also includes the hydrolysis, inert particulate mass balance. The hydrogenotrophic

methanogenesis is explored by including the mass balance of the dissolved hydrogen gas

and the hydrogenotrophic biomass alike. The inclusion of an extensive ionic balance,

featuring ammonia (both ionic and unprotonated form) and inorganic carbon as well as

other cationic and anionic forms serves the purpose to characterize the effluent pH and

corroborate the inorganic nitrogen & carbon balances.

Finally, the inclusion of pH-induced inhibition is useful to describe the possible

process arrest due to the excessive hydrogen or VFA accumulation.

Materials & Methods

43

2.3.1 Biological process overview

The core of the model is comprised by an ODE system of 18 equations. A variable number of

algebraic equations can be present. The table below illustrates the state variables of the

ODE system.

VARIABLE DESCRIPTION MEASURE

X01 Degradable Particulate Matter kgCOD/m3

Xnd Non-degradable Particulate

Matter

kgCOD/m3

X1 Acidogenic Flora kgCOD/m3

X2 Acetoclasic Methanogens kgCOD/m3

X3 Hydrogenotrophic Methanogens kgCOD/m3

S1 Soluble COD kgCOD/m3

S2 VFAS kgCOD/m3

S4 Dissolved Hydrogen kgCOD/m3

S5 Dissolved Methane kgCOD/m3

IN Inorganic Nitrogen M

Snh3 Ammonia M

IC Inorganic Carbon M

Shco3 Hydrogencarbonate M

Scat Cations M

San Anions M

Sgas,h2 Gaseous Hydrogen M

Sgas,co2 Gaseous Carbon Dioxide M

Sgas,ch4 Gaseous Methane M

S3 Dissolved CO2 M

Snh4 Ammonium M

Table 10: State Variables of the modified AMOCO model

Materials & Methods

44

The influent particulate is hydrolyzed to give non-degradable particulate matter (Xnd) and

soluble, non - VFA COD (S1). Such soluble COD is comprised by amino acids, carbohydrate

monomers and long chain fatty acids. This substratum is then processed by the acidogens

(X1) to obtain aqueous hydrogen gas (S4) and VFA (S2): the first is digested by the

hydrogenotrophic methanogens, while the latter by the acetoclastic methanogens to give

dissolved methane (S5). The by-products of the metabolic activities and the bacterial decay

fuel the inorganic nitrogen (IN) and the inorganic carbon (IC) compartments. The acid-base

equilibria-derived mass balances of ammonium and hydrogencarbonate are the

complementary reactions to the inorganic carbon and nitrogen mass. Finally, the equations

of the three components of the gas mix, gaseous carbon dioxide, hydrogen and methane

complete the 18 differential equations. The process diagram (Figure 13) gives a bird’s eye

view over the biological process,

Figure 13: Process Flow Diagram - Biological Process

Materials & Methods

45

The process’ Peterson matrix is a more detailed way to understand the biochemical

dynamic. Its columns display the process’ variables while rows the processes occurring. The

final column displays the rate equations for each process. The matrix is reported in Table 11.

rate

𝐾ℎ1𝑋01

𝜇1 max

𝑌 x1

𝑠 1𝑆 1+𝐾𝑠1I 𝑝𝐻1I 𝐻2X1

Y𝑥2

𝜇2 max

𝑌 x2

𝑠 2𝑆 2+𝐾𝑠2I 𝑝𝐻2X2

Y𝑥3

𝜇3 max

𝑌 x3

𝑠 4𝑆 4+𝐾𝑠42

I 𝑝𝐻3I 𝐼𝐶X3

𝐾d1𝑋1

𝐾d2𝑋2

𝐾d3𝑋3

𝑘𝐿𝑎𝐻2(𝑆 4 16−𝐾𝐻,𝐻2𝑃 𝑔𝑎𝑠,𝐻2)

𝑘𝐿𝑎(𝑆 5 64−𝐾𝐻,𝐶𝐻5𝑃 𝑔𝑎𝑠,𝐶𝐻4)

𝑘𝐿𝑎( 𝑆3−𝐾𝐻,𝐶𝑂2𝑃 𝑔𝑎𝑠,𝐶𝑂2)

𝑘𝐿𝑎𝐻2( 𝐾𝐻,𝐻2𝑃𝑏𝑢𝑏𝑏𝑙𝑒,𝐻2−𝑆 4 64)

∙(−VALVh2)

IN

Nx0

1 -

Ns1

Ns1

-

Y x1N

bac

-Yx2

Nb

ac

-Yx3

Nb

ac

Nb

ac -

Nx0

1

Nb

ac -

Nx0

1

Nb

ac -

Nx0

1

IC

Cx0

1-C

s1

(Cs1

- (Y

x1C

bac

+ (1

-

Yx1

)Cs2

fvfa

)

) (Cs2

-

(Yx2

Cb

ac +

(1

- Y x

2) C

ch4))

(-(Y

x3C

bac

+

(1-Y

x3)

Cch

4))

Cb

ac-C

x01

Cb

ac-C

x01

Cb

ac-C

x01

-1

S 5

1 -

Yx2

1 -

Yx3

-64

S 4

(1 -

Yx1

)fh

2

-1

-16

+16

(1-f

byp

ass)

S 2

(1 -

Yx1

)

f vfa

-1

S 1

1-k

nd

-1

X3

Y x3

-1

X2

Y x2

-1

X1

Y x1

-1

Xn

d

k nd

X0

1

-1

1

1

1

Pro

cess

Dis

inte

grat

ion

an

d

Hyd

roly

sis

S 1 u

pta

ke

S 2 u

pta

ke

S 4 u

pta

ke

X1

dec

ay

X2

dec

ay

X3

dec

ay

H2 g

as

tran

sfer

CH

4 m

ass

tran

sfer

CO

2 m

ass

tran

sfer

Hyd

roge

n

inje

ctio

n

Materials & Methods

46

Table 11: Petersen matrix

2.3.2 Model control

The Model control is achieved by the implementation of a controller that acts upon

determined parameter of the cation, anion, and hydrogen differential equations In this way

desired pH and hydrogen injection levels can be obtained

To regulate these equations, a proportional integrative (PI) controller is implemented. It is a

logical-mathematical operator that, given a target signal �̂�, acts upon the signal that needs

to be controlled, s. If we define an error as (eq. 2.2):

𝑒(𝑡) = �̂�(𝑡) − 𝑠(𝑡) 2.2

the response signal y will be proportional to the sum of the error’s integral over the

simulation time and the error itself. The proportional part on its own makes the controlled

signal to converge to a new value (at steady state conditions) while the integral part ensures

that this value is, in fact the target signal. The constants Kp and KI express the degree and

the rapidity of the response and their value is usually a tradeoff between decreasing

overshoot and increasing convergence time.

𝑦(𝑡) = 𝐾𝑝𝑒(𝑡) + 𝐾𝐼∫ 𝑒(𝜏)𝑑𝜏𝑡

𝑡0

2.3

So, for example, if the target signal is pH = 7.5, the PI controller will open (or close) the

valve parameter of cation/anion differential equations (eqns. 2.20-2.21) with intensity given

by eq. 2.3 till the two signals converge. An example of pH control via PI controller is shown

in the Figure 14, where the PI controller acts upon the anion equation to increase the pH: a

low integral gain results in a high convergence time.

Materials & Methods

47

Figure 14: PI pH Control – the control starts at day 50.

The hydrogen injection is modeled by implementing a negative feedback PI controller that

acts on the VALVH2 parameter of the gas transfer equation 2.4.

𝑟𝑇,𝐻2,𝑏𝑢𝑏𝑏𝑙𝑒 = 𝑘𝑙𝑎,ℎ2 (𝑝ℎ2,𝑏𝑢𝑏𝑏𝑙𝑒𝑘𝐻,ℎ2 −𝑆416)𝑉𝐴𝐿𝑉ℎ2 2.4

The target signal is a succession of increasing steps, as shown in Figure 15.

Figure 15:Hydrogen injection rate as obtained by a PI controller

2.3.3 Equations

The model’s equations are comprised of three main parts: organic matter, the ionic balance

and the gaseous phase.

2.3.3.1 Particulate matter

The equations of the volatile solids are five equations targeting the dynamics of particulate

(both degradable and undegradable), and biomass. The biomass balance is influenced by its

Target pH Controlled

Materials & Methods

48

growth and decay. The growth is modelized with Monod equation, while the decay is a first

order reaction.

𝑑𝑋01𝑑𝑡

=𝑄

𝑉(𝑋01,𝑖𝑛 − 𝑋01) + 𝑘ℎ1𝑋01 + 𝑘𝑑,1𝑋1 + 𝑘𝑑,2𝑋2 + 𝑘𝑑,3𝑋3

2.5

𝑑𝑋𝑛𝑑𝑑𝑡

=𝑄

𝑉(−𝑋𝑛𝑑) + 𝑘𝑛𝑑𝑘ℎ1𝑋01

2.6

𝑑𝑋1𝑑𝑡

=𝑄

𝑉(𝑋1,𝑖𝑛 − 𝑋1) + 𝑌1

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2 − 𝑘𝑑,1𝑋1 2.7

𝑑𝑋2𝑑𝑡

=𝑄

𝑉(𝑋2,𝑖𝑛 − 𝑋2) + 𝑌2

𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2 − 𝑘𝑑,2𝑋2 2.8

𝑑𝑋3𝑑𝑡

=𝑄

𝑉(𝑋3,𝑖𝑛 − 𝑋3) + 𝑌3

𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶 − 𝑘𝑑,3𝑋3 2.9

2.3.3.2 Soluble COD

Soluble COD is comprised by four equations targeting the dynamics of S1 (eq.2.10) VFA (eq

2.11), hydrogen (eq. 2.12) and methane (eq. 2.13). The mass balances are influenced by

biomass’ uptake and gas transfer.

𝑑𝑆1𝑑𝑡=𝑄

𝑉(𝑆1,𝑖𝑛 − 𝑆1) + (1 − 𝑘𝑛𝑑)𝑘ℎ1𝑋01 −

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2 2.10

𝑑𝑆2𝑑𝑡=𝑄

𝑉(𝑆2,𝑖𝑛 − 𝑆2) + (1 − 𝑌𝑥,1)𝑓𝑉𝐹𝐴

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2

−𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2

2.11

𝑑𝑆4𝑑𝑡=𝑄

𝑉(𝑆4,𝑖𝑛 − 𝑆4) + (1 − 𝑌𝑥,1)𝑓𝐻2

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2

−𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶−16𝑘𝑙𝑎,𝐻2 (𝑆416− 𝑘𝐻,ℎ2𝑝ℎ2)

+(1 − 𝑓𝑏𝑦𝑝𝑎𝑠𝑠)16𝑘𝑙𝑎,𝐻2(𝑆416− 𝑘𝐻,ℎ2𝑝𝑏𝑢𝑏𝑏𝑙𝑒)𝑉𝐴𝐿𝑉𝐻2

2.12

𝑑𝑆5𝑑𝑡=𝑄

𝑉(𝑆5,𝑖𝑛 − 𝑆5) + (1 − 𝑌𝑥,1)

𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2

+(1 − 𝑌𝑥,3)𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶

−64𝑘𝑙𝑎 (𝑆564− 𝑘𝐻,𝑐ℎ4𝑝𝑐ℎ4)

2.13

Materials & Methods

49

2.3.3.3 Inorganic carbon

The inorganic carbon equations are two differential and one algebraic equation expressing

the hydrogencarbonate and the carbon dioxide dynamic, The Inorganic carbon dynamic (eq

2.14)) is governed by the CO2 exchange between the biochemical processes while the

hydrogencarbonate equation (eq. 2.15) derives from the acid base reaction equation with

𝐾𝑎𝑏,𝑐𝑜2parameter to account for the fact that the reaction rate is instantaneous. The

algebraic equation (eq. 2.16) expresses the relationship between the two components.

𝑑𝐼𝐶

𝑑𝑡=𝑄

𝑉(𝐼𝐶𝑖𝑛 − 𝐼𝐶) + (𝐶𝑆1 − (1 − 𝑌𝑥,1)𝑓𝑉𝐹𝐴 − 𝑌𝑥,1𝐶𝑏𝑎𝑐)

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2

+(𝐶𝑠2 − (1 − 𝑌𝑥,2)𝐶𝑐ℎ4 − 𝑌𝑥,2𝐶𝑏𝑎𝑐)𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2

+(−(1 − 𝑌𝑥,3)𝐶𝑐ℎ4 − 𝑌𝑥,3𝐶𝑏𝑎𝑐)𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶−𝑘𝑙𝑎(𝑆3 − 𝑘𝐻,𝑐𝑜2𝑝𝑐𝑜2)

+(𝐶𝑥01 − 𝐶𝑥1)𝑘ℎ1𝑋01

2.14

𝑑𝑆ℎ𝑐𝑜3𝑑𝑡

= 𝐾𝑎𝑏,𝑐𝑜2(𝑆ℎ𝑐𝑜3 − (𝑘𝑎,𝑐𝑜2 + 𝑆ℎ) − 𝑘𝑎,𝑐𝑜2𝐼𝐶) 2.15

𝐼𝐶 = 𝑆3 + 𝑆ℎ𝑐𝑜3 2.16

The inorganic Carbon dynamic is simplified because it does not account for all the species.

However, as disputed here, the assumption is not outrageous for a reasonable pH range.

The Inorganic carbon presence in water is caused mainly by the atmospheric CO2 which is

roughly 0.004 % of tropospheric and stratospheric volume. The gaseous carbon dioxide is in

gaseous equilibrium with the dissolved CO2. Such equilibrium is expressed by Henry

constant, which is function of temperature. The aqueous carbon dioxide, on the other hand,

reacts reversibly with water to form carbonic acid. Because carbonic acid is a diprotic acid,

this and will dissociate in hydrogencarbonate (3), and then, to carbonate ion (4).

{

𝐶𝑂2,(𝑔) ↔ 𝐶𝑂2(𝑎𝑞)

𝐶𝑂2(𝑎𝑞) + 𝐻2𝑂(𝑙) ↔ 𝐻2𝐶𝑂3(𝑎𝑞) + 𝐻3𝑂3(𝑎𝑞)+

𝐻2𝐶𝑂3(𝑎𝑞)− + 𝐻2𝑂(𝑙) ↔ 𝐻𝐶𝑂3(𝑎𝑞)

− +𝐻3𝑂3(𝑎𝑞)+

𝐻𝐶𝑂3(𝑎𝑞)− + 𝐻2𝑂(𝑙) ↔ +𝐻3𝑂3(𝑎𝑞)

+ + 𝐶𝑂3(𝑎𝑞)2−

These equilibria lead to the system of equations:

Materials & Methods

50

{

[𝐶𝑂2,(𝑔)] = 𝐾𝑐𝑜2[𝐶𝑂2(𝑎𝑞)]

[𝐻2𝐶𝑂3(𝑎𝑞)][𝐻3𝑂3(𝑎𝑞)+ ]

[𝐶𝑂2(𝑎𝑞)]= 𝐾𝑏

[𝐻𝐶𝑂3(𝑎𝑞)− ][𝐻3𝑂3(𝑎𝑞)

+ ]

𝐻2𝐶𝑂3(𝑎𝑞)− = 𝐾𝑎1

[𝐻3𝑂3(𝑎𝑞)+ ][𝐶𝑂3(𝑎𝑞)

2− ]

[𝐻𝐶𝑂3(𝑎𝑞)− ]

= 𝐾𝑎1

By solving the system

• Expressing the variables as a pH function Kk

• Using the constants at 25 °C and pressure of 1 atm.

• Expressing the variables as a ratio between the species’ concentration and total

inorganic carbon

The following graph is obtained, where the relative abundances of hydrogencarbonate,

carbonate and dissolved carbon dioxide are plotted. The relative abundance of carbonic acid

is also shown (ordinate axis on the left side of the graph).

Figure 16: Inorganic Carbon Speciation

It is evident then, that in the operating pH range the lion’s share of the inorganic carbon is

comprised by the hydrogencarbonate and dissolved carbon dioxide. It is then assumed that

these are the only inorganic species present, thus simplifying the problem.

Materials & Methods

51

2.3.3.4 Inorganic nitrogen

Inorganic nitrogen is the sum of NH3 and NH4+. The dynamic of these species is influenced by

the influent and the biomass activity, while the speciation between the protonated and

unprotonated form is influenced by pH and temperature.

Figure 17: Un-ionized ammonia percentage as a function of pH and temperature (Speece, 1996)

Ammonium cation is used as a macronutrient source for bacterial growth. It is a common

practice to assume that the brute biomass’ formula is C5H7O2N; it follows then that the

nitrogen is on the average 2 % of the cell’s mass composition, being absorbed during the

growth and released in the decay process. The main reasons to include the nitrogen

equations is to have a grasp on the nitrogen balance and to foresee possible ammonia-

induced inhibition and to evidence possible events of nitrogen-induced starvation.

The inorganic nitrogen equation eq. 2.17 is obtained either from the input, the biomass

decay or the processing of digestible particulate and soluble COD. The relation between the

ammonia and the ammonium cation is expressed by (eq. 2.18) where 𝐾𝑎𝑏,𝐼𝑁 is an expedient

to transform an algebraic expression into a differential one. The nitrogen mass balance is

expressed by the algebraic equation (eq. 1.19)

Materials & Methods

52

𝑑𝐼𝑁

𝑑𝑡=𝑄

𝑉(𝐼𝑁𝑖𝑛 − 𝐼𝑁) + (𝑁𝑥01 − 𝑁𝑥1)𝑘ℎ1𝑋01 + (𝑁𝑠1 − 𝑌𝑥1𝑁𝑏𝑎𝑐)

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2

−𝑌𝑥,2𝑁𝑏𝑎𝑐𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2

−𝑌𝑥,3𝑁𝑏𝑎𝑐𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶 + (𝑁𝑏𝑎𝑐 + 𝑁𝑥01)(𝑘𝑑,1𝑋1 + 𝑘𝑑,2𝑋2 + 𝑘𝑑,3𝑋3)

2.17

𝑑𝑆𝑛ℎ3𝑑𝑡

= 𝐾𝑎𝑏,𝐼𝑁(𝑆𝑛ℎ3 − (𝑘𝑎,𝐼𝑁 + 𝑆ℎ) − 𝑘𝑎,𝐼𝑁𝐼𝑁)

2.18

𝐼𝑁 = 𝑆𝑛ℎ3 + 𝑆𝑛ℎ4 2.19

2.3.3.5 Anions and cations

The anionic and cationic equations are featured to model a possible inflow and intestine

development of strong bases and acids and are useful to control the model’s pH.

𝑑𝑆𝑐𝑎𝑡𝑑𝑡

=𝑄

𝑉(𝑆𝑐𝑎𝑡,𝑖𝑛 − 𝑆𝑐𝑎𝑡) +

𝑞𝑎𝐶𝑎𝑉𝑉𝐴𝐿𝑉𝑎 2.20

𝑑𝑆𝑎𝑛𝑑𝑡

=𝑄

𝑉(𝑆𝑎𝑛,𝑖𝑛 − 𝑆𝑎𝑛) +

𝑞𝑏𝐶𝑏𝑉𝑉𝐴𝐿𝑉𝑏 2.21

The model’s pH is regulated by the VALV parameters which can be user defined or

mimicking a data series.

2.3.3.6 Gas equations

The gas equations comprised by a set of three differential equations and a variable number

of ancillary algebraic equations concerning the gas outflow, transfer rates and pressures.

The methane (eq. 2.23) and carbon dioxide (eq. 2.24) balances are influenced solely by the

transfer from the species’ aqueous phase, while the hydrogen (eq. 2.22) also encompasses a

flow derived from bypass of the bubbling process. The equation 2.25 expresses the gas

outflow, which is influenced solely by the gas transfer from liquid phase.

𝑑𝑆𝑔𝑎𝑠,ℎ2

𝑑𝑡= −𝑞𝑔𝑎𝑠

𝑆𝑔𝑎𝑠,ℎ2

𝑉𝑔𝑎𝑠+𝑉

𝑉𝑔𝑎𝑠𝑘𝑙𝑎,ℎ2(

𝑆416− 𝑘𝐻,ℎ2𝑝ℎ2)

+ 𝑓𝑏𝑦𝑝𝑎𝑠𝑠𝑘𝑙𝑎,ℎ2(𝑘𝐻,ℎ2𝑝ℎ2,𝑏𝑢𝑏𝑏𝑙𝑒 −𝑆416)

2.22

Materials & Methods

53

𝑑𝑆𝑔𝑎𝑠,𝑐ℎ4

𝑑𝑡= −𝑞𝑔𝑎𝑠

𝑆𝑔𝑎𝑠,𝑐ℎ4

𝑉𝑔𝑎𝑠+𝑉

𝑉𝑔𝑎𝑠𝑘𝑙𝑎 (

𝑆564− 𝑘𝐻,𝑐ℎ4𝑝𝑐ℎ4) 2.23

𝑑𝑆𝑔𝑎𝑠,𝑐𝑜2

𝑑𝑡= −𝑞𝑔𝑎𝑠

𝑆𝑔𝑎𝑠,𝑐𝑜2

𝑉𝑔𝑎𝑠+𝑉

𝑉𝑔𝑎𝑠𝑘𝑙𝑎(𝑆3 − 𝑘𝐻,𝑐𝑜2𝑝𝑐𝑜2)

2.24

𝑞𝑔𝑎𝑠 =𝑅𝑇

(𝑝𝑔𝑎𝑠 − 𝑝𝑣)𝑉(𝑘𝑙𝑎,ℎ2(

𝑆416− 𝑘𝐻,ℎ2𝑝ℎ2) + 𝑘𝑙𝑎 (

𝑆564− 𝑘𝐻,𝑐ℎ4𝑝𝑐ℎ4) + 𝑘𝑙𝑎(𝑆3

− 𝑘𝐻,𝑐𝑜2𝑝𝑐𝑜2))

2.25

2.3.3.7 pH equations

The pH is expressed by an algebraic equation which derives from the charge balance, by

expressing the hydroxide concentration in terms of hydronium and the solving the

substitution variable 𝜃 the resulting quadratic expression.

𝜃 = 𝑆𝑐𝑎𝑡 + 𝑆𝑛ℎ4 − 𝑆ℎ𝑐𝑜3 −𝑆216− 𝑆𝑎𝑛 2.26

𝑆ℎ =−𝜃 + √𝜃2 − 4𝑘𝑤

2

2.27

𝑝𝐻 = −𝑙𝑜𝑔(𝑆ℎ) 2.28

2.3.3.8 Inhibition equations

The inhibition is either pH, H2 or inorganic carbon induced. The deleterious effects of acidic

environments upon the microbial growth are described by the pH inhibition equations (eq.

2.29-2.31), while the inhibition of hydrogenotrophes by low levels of inorganic carbon is

described by eq. 2.33. The acidogenic fauna on the other hand is negatively influenced by

high levels of hydrogen (eq. 2.32). The whole inhibition dynamic is useful to describe the

process collapse in the case of excessive hydrogen injection leading on one hand to

excessive alkalinity consumption and subsequent hydrogenotrophic starvation, while on the

other the collapse of acidogens, ill-tolerating high hydrogen concentration.

IpH1 = {exp(−3(

pH − pHUL,1pHUL,1 − pHLL,1

)

2

) ∶ pH < pHUL.1

1 ∶ pH > pHUL,1

2.29

Materials & Methods

54

IpH2 = {exp(−3(

pH − pHUL,2rHUL,2 − pHLL,2

)

2

) ∶ pH < pHUL,2

1 ∶ pH > pHUL,2

2.30

lpH3 = {exp(−3(

pH − pHUL,3pHUL,3 − pHLL.3

)

2

) ∶ pH < pHUL,3

1 ∶ pH > pHUL,3

2.31

IH2 =1

1 +S4KI,H2

2.32

𝐼𝐼𝐶 =1

1 +𝐾𝐼,𝐼𝐶𝐼𝑐

2.33

The pH inhibition function is a composition between a gaussian and a constant function.

Observation must be drawn to the fact that the pH-induced inhibition is active only for the

acidic environments, so the reader must employ a careful approach when modeling at high

pH. It is advised to maintain the pH in 6.5-8.2 range (Speece, 1996).

Figure 18: IpH function illustration - comparison with a pure gaussian function

2.3.4 Parameters

The parameters of the model are obtained mainly from the works of (Rosen &

Jeppesen 2006) but some of the parameters were obtained by (Rosito, 2019): namely the

hydrolysis rate and the particulate undegradable fractions. Most of the input parameters

-0.2

0

0.2

0.4

0.6

0.8

1

1.2

0 2 4 6 8 10 12 14 16

pH

pure gaussian composite gaussian

Materials & Methods

55

were calculated from the experiment’s data, while inorganic nitrogen and soluble COD

inputs were calculated from values obtained from Bresso plant sludge analysis (Table 12).

Parameter Value Unit

Alkalinity 1500 𝑚𝑔𝐶𝑎𝐶𝑂3

/𝑙 Total particulate COD 18 𝑔𝐶𝑂𝐷 𝑙⁄

pH 6.67

VFAs (total eq. of acetate) 0.48078 𝑔𝐶𝑂𝐷 𝑙⁄

Ammonium ion 253

mg/L

𝑚𝑔𝑁𝐻4+ 𝑙⁄

Total soluble COD

Sol

4.7 𝑔𝐶𝑂𝐷 𝑙⁄ Table 12: Bresso sludge characterization (Rosito, 2018) -

2.3.4.1 Input parameters

Since the data does not cover all the model’s input species, some of the parameters must be

calculated.

Substrate & biomass

The data (Table 12) for the input volatile solids is available. It is a time-dependent

signal because substrate comes from multiple tanks that are changed throughout the

experiment. The variability of this signal can be observed in Figure 19.

Figure 19: particulate COD input

The input sludge is therefore a time dependent function.

𝑋01,𝑖𝑛 = 𝑓(𝑡) [𝑘𝑔𝐶𝑂𝐷

𝑚3]

2.34

10

15

20

25

30

35

40

0 20 40 60 80 100 120 140 160

particulate COD inputgCOD/l

Materials & Methods

56

In regard to the soluble substrate, since there was no available data, it was decided to

express their concentration in proportion to particulate COD. By using the data in Table 12,

the following ratios were defined for the total soluble COD and the VFA:

𝑅𝐶𝑂𝐷,𝑠𝑜𝑙 =𝐶𝑂𝐷𝑠𝑜𝑙,𝑡𝑜𝑡𝐶𝑂𝐷𝑋,𝑇𝑂𝑇

2.35

𝑅𝐶𝑂𝐷,𝑉𝐹𝐴 =𝐶𝑂𝐷𝑠𝑜𝑙,𝑉𝐹𝐴𝐶𝑂𝐷𝑋,𝑇𝑂𝑇

2.36

On the basis of these ratios, new values were calculated:

𝑆2,𝑖𝑛 = 𝑅𝐶𝑂𝐷,𝑉𝐹𝐴𝑋01,𝑖𝑛 2.37

𝐶𝑂𝐷𝑆𝑜𝑙,𝑇𝑂𝑇,𝑖𝑛 = 𝑅𝐶𝑂𝐷,𝑠𝑜𝑙𝑋01,𝑖𝑛 2.38

𝑆1,𝑖𝑛 = 𝐶𝑂𝐷𝑆𝑜𝑙,𝑇𝑂𝑇,𝑖𝑛 − 𝑆2,𝑖𝑛 2.39

In regard to the biomass input, were assumed the same values as hypothesized by (Rosen &

Jeppesen, 2006)

Inorganic Carbon

The inorganic carbon concentration is calculated by using the alkalinity data. The

alkalinity is a way to measure the solution’s ability to buffer acid addition. It is therefore

expressed by the sum of the concentrations of ionic species that are able to react with

hydronium. The equation 2.40 can be solved to obtain hydrogencarbonate concentration.

𝐴𝑙𝑘 = 𝑆2,𝑖𝑜𝑛 + 𝑆ℎ𝑐𝑜3− + 𝑆𝑂𝐻 + 𝑆𝑛ℎ3 − 𝑆ℎ 2.40

Where the hydroxide and acetate concentrations are expressed as:

Materials & Methods

57

𝑆2,𝑖𝑜𝑛 =

𝑆264 𝑘𝑎,𝑎𝑐

𝑘𝑎,𝑎𝑐 + 𝑆ℎ

2.41

𝑆𝑜ℎ =𝑘𝑤𝑆ℎ

2.42

The inorganic carbon is therefore obtained from the hydrogencarbonate – carbon dioxide

equilibrium:

𝐶𝑂2(𝑎𝑞) + 𝐻2𝑂(𝑙)𝑘𝑎,𝑎𝑐↔ 𝐻𝐶𝑂3(𝑎𝑞)

− + 𝐻3𝑂(𝑎𝑞)+

2.43

{

𝐼𝐶 = 𝑆3 + 𝑆ℎ𝑐𝑜3−

𝐼𝐶 =𝑆ℎ𝑐𝑜3𝑘𝑎,𝑐𝑜2𝑘𝑎,𝑐𝑜2 + 𝑆ℎ

2.44

Inorganic Nitrogen

The Inorganic nitrogen input was calculated on the basis of the Brosso sludge data (Table

12), where ammonium cation molarity is known. It was assumed that throughout all the

experiment’s duration the input ammonium was constant. To determine the molarity of

ammonia, the acid base dissociation equilibrium as bell as mass balance equations were

used. By writing the acid-base reaction of ammonia and by seeing the ammonium as a weak

acid we obtain:

𝑁𝐻4(𝑎𝑞)+ + 𝐻2𝑂(𝑙)

𝑘𝑎,𝐼𝑁↔ 𝑁𝐻3(𝑎𝑞) + 𝐻3𝑂(𝑎𝑞)

+ 2.45

And by applying the inorganic nitrogen mass balance following system is obtained:

Materials & Methods

58

{

𝑆𝑛ℎ3𝑆ℎ𝑆𝑛ℎ4

= 𝑘𝑎,𝐼𝑁

𝐼𝑁 = 𝑆𝑛ℎ3 + 𝑆𝑛ℎ4

2.46

Thus, the inorganic nitrogen calculated:

𝐼𝑁 =𝑆𝑛ℎ4(𝑘𝑎,𝐼𝑁 + 𝑆ℎ)

𝑆ℎ

2.47

Anions & cations

The cation input value vas assumed the same as (Rosen & Jeppsen, 2006) while the anion

value was calculated from the ionic balance:

𝑆ℎ + 𝑆𝑛ℎ4 + 𝑆𝑐𝑎𝑡 − (𝑆𝑎𝑛 + 𝑆𝑂𝐻 + 𝑆ℎ𝑐𝑜3 + 𝑆2,𝑖𝑜𝑛) = 0 2.48

Obtained values

The input parameters are therefore calculated and reported in Table 13.

Parameter Value Unit

𝑋01,𝑖𝑛 f(t) 𝑘𝑔𝐶𝑂𝐷 𝑚3⁄

𝑋1,𝑖𝑛 0.01 𝑘𝑔𝐶𝑂𝐷 𝑚3⁄

𝑋2,𝑖𝑛 0.01 𝑘𝑔𝐶𝑂𝐷 𝑚3⁄

𝑋3,𝑖𝑛 0.01 𝑘𝑔𝐶𝑂𝐷 𝑚3⁄

𝑆1,𝑖𝑛 g(t) 𝑘𝑔𝐶𝑂𝐷 𝑚3⁄

𝑆2,𝑖𝑛 h(t) 𝑘𝑔𝐶𝑂𝐷 𝑚3⁄

𝐼𝐶𝑖𝑛 0.006979 𝑀

𝐼𝑁𝑖𝑛 0.014108 𝑀

𝑆𝑐𝑎𝑡,𝑖𝑛 0.02 𝑀

𝑆𝑎𝑛,𝑖𝑛 0.018127 𝑀

Table 13: Input parameters

Materials & Methods

59

2.3.4.2 Biochemical parameters

The biochemical parameters deal with the description of the biomass growth, decay, and

hydrolysis.

Parameter Value Unit Description

𝜇1,𝑚𝑎𝑥 0.3 1 𝑑⁄ Max growth rate - acidogens

𝜇2,𝑚𝑎𝑥 0.33 1 𝑑⁄ Max growth rate - acetoclasts

𝜇3,𝑚𝑎𝑥 2 1 𝑑⁄ Max growth rate - hydrogenotrophes

𝑘𝑑,1 0.02 1 𝑑⁄ Decay rate - acidogens

𝑘𝑑,2 0.02 1 𝑑⁄ Decay rate - acetoclasts

𝑘𝑑,3 0.02 1 𝑑⁄ Decay rate – hydrogenotrophes

𝐾𝑠,1 0.5 𝑘𝑔𝐶𝑂𝐷𝑚3

Half saturation constant for S1 uptake

𝐾𝑠,2 0.04 𝑘𝑔𝐶𝑂𝐷𝑚3

Half saturation constant for VFA uptake

𝐾𝑠,4 7e-6 𝑘𝑔𝐶𝑂𝐷𝑚3

Half saturation constant for hydrogen

uptake

𝑌𝑥,1 0.15 𝑘𝑔𝐶𝑂𝐷,𝑥𝑘𝑔𝐶𝑂𝐷,𝑠

Yield coefficient – acidogens

𝑌𝑥,2 0.05 𝑘𝑔𝐶𝑂𝐷,𝑥𝑘𝑔𝐶𝑂𝐷,𝑠

Yield coefficient – acetoclasts

𝑌𝑥,4 0.06 𝑘𝑔𝐶𝑂𝐷,𝑥𝑘𝑔𝐶𝑂𝐷,𝑠

Yield coefficient – hydrogenotrophes

𝑘ℎ1 0.3855 1 𝑑⁄ Hydrolysis rate

𝑘𝑛𝑑 0.599 𝑘𝑔𝐶𝑂𝐷,𝑥𝑘𝑔𝐶𝑂𝐷,𝑥

Undegradable fraction of volatile solids

Table 14: Biochemical parameters

Materials & Methods

60

2.3.4.3 Physiochemical parameters

These parameters describe the physiochemical properties of substances such as acid-base

reactions, gas mass transfer and water self-ionization.

Parameter Value Unit Description

𝑘𝑙𝑎 8.3 1 ℎ⁄

Methane and carbon

dioxide gas transfer

rate

𝑘𝑙𝑎,ℎ2 1..75 1 ℎ⁄ Hydrogen gas

transfer rate

𝑝𝑣 0.0313 ∙ 𝑒𝑥𝑝 [[5290 ((1

298−1

𝑇))]] 105 𝑃𝑎 Vapor pressure

𝑘𝐻,ℎ2 0.035 ∙ 𝑒𝑥𝑝 [−19410

𝑅(1

298−1

𝑇)] 105

𝑀

𝑃𝑎

Hydrogen Henry

constant

𝑘𝐻,𝑐ℎ4 0.0014 ∙ 𝑒𝑥𝑝 [−14290

𝑅(1

298−1

𝑇)] 105

𝑀

𝑃𝑎

Methane henry

constant

𝑘𝐻,𝑐𝑜2 0.00076 ∙ 𝑒𝑥𝑝 [−4180

𝑅(1

298−1

𝑇)] 105

𝑀

𝑃𝑎

Carbon dioxide henry

constant

𝑘𝑤 𝑒𝑥𝑝 [148.98 − 13847.26

𝑇 − 23.65 ∗ 𝑙𝑜𝑔(𝑇)] 𝑀2

Water self-ionization

constant

𝐾𝑎𝑏,𝑐𝑜2 1010 1

𝑀 ∙ 𝑑

Inorganic carbon

dissociation rate

𝐾𝑎𝑏,𝐼𝑁 1010 1

𝑀 ∙ 𝑑

Inorganic nitrogen

dissociation rate

𝑘𝑎,𝑐𝑜2 10−6.35𝑒𝑥𝑝 [7646

𝑅(1

298−1

𝑇)] 𝑀

Inorganic carbon acid

dissociation constant

𝑘𝑎,𝐼𝑁 10−9.25𝑒𝑥𝑝 [51965

𝑅(1

298−1

𝑇)] 𝑀

Inorganic nitrogen

acid dissociation

constant

𝑘𝑎,𝑎𝑐 10−4.86 𝑀 Acetic acid

dissociation constant

Table 15: Physiochemical parameters

Materials & Methods

61

2.3.4.4 Reactor parameters

These parameters describe the reactor’s operational conditions

Parameter Value Unit Description

Q 0.0005 𝑚3 𝑑⁄ Sludge flow

V 11e-3 m3 Reactor volume

Vgas 3e-3 m3 Reactor headspace volume

pgas 101325 Pa Reactor’s operational pressure

T 315.15 K Temperature

Table 16: Reactor parameters

2.3.4.5 Fractioning parameters

These parameters describe the nitrogen and carbon content of the sludge constituents as

well as the soluble COD fractions that become hydrogen and VFA. The injected hydrogen

bypass parameter can also be found here.

Parameter Value Unit Description

𝐶𝑥01 0.02786 𝑘𝑚𝑜𝑙𝐶𝑚3

Volatile solids carbon content

𝐶𝑠1 0.0313 𝑘𝑚𝑜𝑙𝐶𝑚3

Soluble COD carbon content

𝐶𝑠2 0.0313 𝑘𝑚𝑜𝑙𝐶𝑚3

VFA carbon content

𝐶𝑐ℎ4 0.0156 𝑘𝑚𝑜𝑙𝐶𝑚3

Methane carbon content

𝐶𝑏𝑎𝑐 0.0313 𝑘𝑚𝑜𝑙𝐶𝑚3

Biomass’ carbon content

𝑁𝑥01 0.002 𝑘𝑚𝑜𝑙𝑁𝑚3

Volatile solids nitrogen

content

𝑁𝑠1 0.007 *

0.3

𝑘𝑚𝑜𝑙𝑁𝑚3

Soluble COD nitrogen content

𝑁𝑏𝑎𝑐 0.00625 𝑘𝑚𝑜𝑙𝑁𝑚3

Biomass’ nitrogen content

𝑓𝑣𝑓𝑎 0.482254 [/] Fraction of S1 that is processed

to VFA

Materials & Methods

62

𝑓ℎ2 0.517746 [/] Fraction of S1 that is processed

to H2

𝑓𝑏𝑦𝑝𝑎𝑠𝑠 0.05 [/] Injected hydrogen bypass to

headspace

Table 17: Fractioning parameters

2.3.4.6 Inhibition parameters

These parameters describe the inhibition functions. The pH inhibition parameters describe

the conditions below which maximum inhibition occurs (pHll) and above which there is no

inhibition (pHul).

Parameter Value Unit Description

𝑝𝐻𝑢𝑙,1 4.5 [/] X1 growth inhibition upper limit

𝑝𝐻𝑙𝑙,1 5.5 [/] X1 growth inhibition lower limit

𝑝𝐻𝑢𝑙,2 5.8 [/] X2 growth inhibition upper limit

𝑝𝐻𝑙𝑙,2 6.7 [/] X2 growth inhibition lower limit

𝑝𝐻𝑢𝑙,3 5 [/] X3 growth inhibition upper limit

𝑝𝐻𝑙𝑙,3 6 [/] X3 growth inhibition lower limit

𝐾𝐼,𝐼𝐶 1e-4 𝑀

Inorganic carbon inhibition

constant

𝐾𝐼,ℎ2 1e-5 𝑘𝑔𝐶𝑂𝐷𝑚3

Hydrogen inhibition constant

Table 18: Inhibition parameters

Materials & Methods

63

2.4 LFT MODEL

The parametric identification of the system is complex because the system is nonlinear. To

overcome this problematic an algorithm was developed: the Linear Fractional Transformation.

The Linear Fractional Transformation (LFT) is a methodology of parametric identification for

nonlinear systems that was developed by (Ferretti & Della Bona. 2015) The main idea behind

this method is the linearization of the nonlinear parts by the means of substitution of non-

linear snippets of the equation system with variables. For each iteration the array of non-

linear variables is calculated and the then the linearized system is solved with common

optimization algorithms.

Given a non-linear, time-invariant, multi input multi-output (MIMO), time-continuous

system:

where x ϵ ℝn, y ϵ ℝn, are the state and noise-free output vectors, with u ϵ ℝm being the

input vector,

and 𝛿° = {𝛿°𝑖} (i=1,…q} a vector of unknown parameters. The system can be formulated in

LFT form:

where z ϵ ℝn, 𝜔 ϵ ℝn , w ϵ ℝn , 𝜉 ϵ ℝn are vectors of auxiliary variables, A, Bi, Ci, Dij are 16

known constant matrices, ri are the sizes of the corresponding identity matrices 𝐼𝑟𝑖, in the Δ

block and 𝛩(𝜔(𝑡)):ℝ𝑛→ ℝ𝑛 is a known nonlinear vector function.

To aid the derivation of such formulation, a symbolic computing approach was developed.

The model, written in the LFT formulation can be divided into three parts:

Materials & Methods

64

• A linear part

• Non-linear part: the theta matrix, which contains snippets of non- linear parts

• Uncertain part: the equation of w(t) where z(t) multiplies the matrix of the unknown

parameters.

The rewriting of the original model allows a direct computation of the Hessian and gradient

of the cost function considered for the identification.

Figure 8: A scheme of the LFT Formulation

The output observation equation is assumed to take the form:

�̂�(𝑡𝑘) = 𝑦(𝑡𝑘) + 휀(𝑡𝑘)

Where 𝑡𝑘, k=1,…..,N denotes the sampling instant and 휀1(𝑡𝑘) is a discrete-time white noise

of the variance σi2. Then, the identification problem can be formulated as follows: given the

sampled data {u(𝑡𝑘), �̌�(𝑡𝑘)} 𝐾=1 find the values of parameters �̌� minimizing the cost function.

Where e(𝑡𝑘,𝛿°)=�̌�(𝑡𝑘)− �̌�(𝑡𝑘,𝛿 ) is the prediction error between the measured output �̌�(𝑡𝑘)

and the output �̌�(𝑡𝑘,𝛿 ) predicted by model, with parameters 𝛿 instead of the true

parameters 𝛿° and W is a weight matrix.

Materials & Methods

65

It is also worth mentioning, that, in order to deal with parameter estimation, it is essential

to rewrite model by introducing normalized unknown parameters 𝛿𝑖 ̅varying between ±1 as

parameters 𝛿𝑖 varies between a maximum 𝛿𝑖,𝑚𝑎𝑥 and a minimum 𝛿𝑖,𝑚𝑖𝑛 with:

The δ̂̌ is a maximum-likehood estimate of the model parameters 𝜹 for output-error and can

be obtained through well-known optimization procedures such as the Gauss-Newton

algorithm

Where v is the iteration number, 𝛼(𝑣) is the step size, 𝑔(�̌� (𝑣)): ℝ𝑛→ ℝ𝑛 and 𝑔(�̌� (𝑣)):ℝ𝑞→ ℝ𝑞

and �̂�(�̌� (𝑣)):ℝ𝑞→ ℝ𝑞𝑥𝑞 are respectively the gradient vector and a positive semi-definite

approximation of the Hessian of the cost function with respect to the unknown parameters:

where 𝐸(𝑘,𝛿)∈ℝ𝑝𝑥𝑞 is the Jacobian of 𝑒(𝑘,𝛿) and is given by:

In this parameter identification, the optimization is carried out from MATLAB fmincon

function, which receives as inputs parameters the cost function 𝐽(𝛿), the gradient 𝑔(𝜹), the

approximated Hessian �̂�(𝜹) and the limits of parameters as restrictions of the optimization

problem.

Materials & Methods

66

2.5 AMOCO MODEL LINEARIZATION

To obtain the AMOCO model in a linearized form a MATLAB script was developed. The aim

of the script is the creation of the LFT matrices starting from the Source Code as written in

Modelica. The script, as reported in the Appendix (link) uses the MATLAB symbolic package

operate the transformation. At first the arrays of state variables, inputs and unknown

parameters are defined

𝒙 =

(

𝑥1𝑥2𝑥3𝑥4𝑥5𝑥6𝑥7𝑥8𝑥9𝑥10𝑥11𝑥12𝑥13𝑥14𝑥15𝑥16𝑥17𝑥18)

=

(

𝑋01𝑋𝑛𝑑𝑋1𝑋2𝑋3𝑆1𝑆2𝑆4𝑆5𝐼𝐶𝐼𝑁𝑆ℎ𝑐𝑜3𝑆𝑛ℎ3𝑆𝑐𝑎𝑡𝑆𝑎𝑛𝑆𝑔𝑎𝑠,ℎ2𝑆𝑔𝑎𝑠,𝑐ℎ4𝑆𝑔𝑎𝑠,𝑐𝑜2)

𝒖 =

(

𝑢1𝑢2𝑢3𝑢4𝑢5𝑢6𝑢7𝑢8𝑢9𝑢10)

=

(

𝑋01,𝑖𝑛𝑋1,𝑖𝑛𝑋2,𝑖𝑛𝑋3,𝑖𝑛𝑆1,𝑖𝑛𝑆2,𝑖𝑛𝐼𝐶𝑖𝑛𝐼𝑁𝑖𝑛𝑆𝑐𝑎𝑡,𝑖𝑛𝑆𝑎𝑛,𝑖𝑛)

𝜹 = (

𝛿1𝛿2𝛿3

) = (

𝑓𝑉𝐹𝐴𝑘𝑙𝑎𝑘𝑙𝑎,ℎ2

)

The array of the non- linear snippets is then written, the ξ or XI array.

𝝃 =

(

𝜉1𝜉2𝜉3𝜉4𝜉5𝜉6𝜉7𝜉8𝜉9𝜉10𝜉11𝜉12𝜉13𝜉14𝜉15𝜉16)

=

(

ω3ω6K𝑠,1 + ω6ω4ω7

K𝑠,2 + ω7ω5ω8

K𝑠,4 + ω8

ω12(−1

2(ω7

64− ω11 + ω12 + ω13 − ω14 + ω15) +

1

2√(

ω7

64− ω11 + ω12 + ω13 − ω14 + ω15)

2− 4𝑘𝑤)

ω13(−1

2(ω7

64− ω11 + ω12 + ω13 − ω14 + ω15) +

1

2√(

ω7

64− ω11 + ω12 + ω13 − ω14 + ω15)

2− 4𝑘𝑤)

ω16(ω816− 𝑅 𝑇 𝐾𝐻,16,ω16)

ω17(ω816− 𝑅 𝑇 𝐾𝐻,16,ω16)

ω18(ω816− 𝑅 𝑇 𝐾𝐻,16,ω16)

ω16(ω964− 𝑅 𝑇 𝐾𝐻,𝐶𝐻4,ω17)

ω17(ω964

− 𝑅 𝑇 𝐾𝐻,𝐶𝐻4,ω17)

ω18(ω964

− 𝑅 𝑇 𝐾𝐻,𝐶𝐻4,ω17)

−ω16(ω12 − ω10 − 𝑅 𝑇 𝐾𝐻,𝑐𝑜2,ω18)

−ω17(ω12 − ω10 − 𝑅 𝑇 𝐾𝐻,𝑐𝑜2,ω18)

−ω18(ω12 − ω10 − 𝑅 𝑇 𝐾𝐻,𝑐𝑜2,ω18)

(−1

2(ω7

64− ω11 + ω12 + ω13 − ω14 + ω15) +

1

2√(

ω7

64− ω11 + ω12 + ω13 − ω14 + ω15)

2− 4𝑘𝑤)

1

(−1

2(ω7

64− ω11 + ω12 + ω13 − ω14 + ω15) +

1

2√(

ω7

64− ω11 + ω12 + ω13 − ω14 + ω15)

2− 4𝑘𝑤) )

Materials & Methods

67

These snippets are then substituted in the non- linear system. The final step is the creation

of the w array, which contains all the variables that are multiplied by the unknown

parameters. This step is mediated by the vector of the z

𝒛 =

(

𝑧1𝑧2𝑧3𝑧4𝑧5𝑧6𝑧7𝑧8𝑧9𝑧10𝑧11𝑧12𝑧13)

=

(

𝜉1𝑥816− RT 𝐾𝐻,ℎ2𝑥16

𝜉6𝜉7𝜉8

𝑥964− RT 𝐾𝐻,𝑐ℎ4𝑥17

𝑥10 − 𝑥12 − RT 𝐾𝐻,𝑐𝑜2𝑥18𝜉9𝜉10𝜉11𝜉12𝜉13𝜉14

)

𝒘 =

(

𝑤1𝑤2𝑤3𝑤4𝑤5𝑤6𝑤7𝑤8𝑤9𝑤10𝑤11𝑤12𝑤13)

=

(

𝛿1𝑧1𝛿2𝑧2𝛿2𝑧3𝛿2𝑧4𝛿2𝑧5𝛿3𝑧6𝛿3𝑧7𝛿3𝑧8𝛿3𝑧9𝛿3𝑧10𝛿3𝑧11𝛿3𝑧12𝛿3𝑧13)

At last, the array of the output variables y is built.

𝒚 =

(

𝑦1𝑦2𝑦3𝑦4𝑦5𝑦6𝑦7𝑦8𝑦9)

=

(

𝑋𝑇𝑂𝑇𝑆1𝑆2𝐴𝑙𝑘𝑞𝑔𝑎𝑠,𝑆𝑔𝑎𝑠,ℎ2𝑆𝑔𝑎𝑠,𝑐ℎ4𝑆𝑔𝑎𝑠,𝑐𝑜2𝑆ℎ )

=

(

𝑥1+𝑥2 + 𝑥3+𝑥4+𝑥5𝑥6𝑥7

(𝑥764 − 𝜉15 + 𝑥10 + 𝑥11 − 𝑥13 + 𝑘𝑤𝜉16)50000

103𝑅𝑇𝑝𝑔𝑎𝑠−𝑝𝐻2𝑂

𝑉(𝑤2 + 𝑤6 + 𝑤7)

𝑥16𝑥17𝑥18𝜉15

)

Chapter 3

3 RESULTS AND DISCUSSION

In this chapter the results of the AMOCO model validation and LFT parametric

identification are discussed. The aim of the model’s validation is to evidence behavior

anomalies and evaluate the degree of similarity between the model and the data. At first a

preliminary system response analysis was performed which aimed to check the correct

functioning and underline model’s possible critical issue. This was obtained by modulating

the substrate input signals with a succession of increasing step signals and then checking the

mass balances of COD, nitrogen, and carbon. This has evidenced several issues such as an

excessive soluble substrate consumption,

To further investigate these issues and to simulate a scenario of VFA accumulation

the parametric set of the biochemical dynamic of the acetoclastic methanogens was

manipulated in such a way as to get VFA values more consistent with the data. Then the

stress response to growing levels of hydrogen injection were evaluated once again

underlining the critical issues of excessive hydrogen consumption. With these issues in

mind, finally the comparison between the model and the available data was made to truly

understand if the model is an adequate predictor of the real behaviour.

Results and Discussion

69

3.1 SLUDGE FLOW VARIATION

To have a preliminary understanding of the model’s behavior the system was

stimulated with the feed overflow and then the response was evaluated and then the mass

balances of COD, nitrogen and carbon were performed. The feed overflow was simulated by

modulating the input of the sludge as a succession of ever-increasing step signals for the

three main components of the inflow sludge.

Figure 20: Sludge Input Signal

time S1,in S2,in X01,in

0 6.27 0.48 18

40 12.65 1.44 54

80 18.91 2.16 81

[d] [gCOD/l] [gCOD/l] [gCOD/l]

Table 19: Concentrations of the Input Species

Results and Discussion

70

3.1.1 Total Solid Output

The volatile solids’ response is significant especially with regards to inert solids. This is

consistent with aim of the sludge stabilization, as shown in the figure below.

Figure 21: Solid Matter Output

3.1.2 Soluble COD output

The soluble COD modulation showed an important problematic within the model.

Regardless of the input concentration the outflow of the soluble COD shows always

extremely low values. As it will be explained in the next sections, this is inconsistent with the

real data. The Figure 22 shows the response signal with comparison to the initial conditions,

while the Figure 23 gives a more precise idea of how weak the output signal dynamic is.

Figure 22: Soluble COD output – comparison with the initial conditions

Results and Discussion

71

Figure 23: Soluble COD output - up close

3.1.3 Biomass

The output signal of acidogens and acetogens (X1) is particularly responsive to the

substrate increase, mainly due to a high Yield Coefficient, an order of magnitude greater

than that of the other groups (see chapter 2). More than that, the model does not consider

possible substratum-induced inhibitions.

Figure 24: Biomass Output Concentration

3.1.4 Inorganic Carbon

The increased biomass activity leads to a consistent accumulation of Inorganic Carbon,

which in turn leads to an increase of pH, as shown in Figure 25 & Figure 26. The signal,

Results and Discussion

72

starting at steady state condition of 6.9, reaches the final value of 7.5. This result also

showcases a weak influence of the substratum influent over the pH, which increases only

slightly, from 6.9 to 7.4.

Figure 25: Inorganic Carbon

3.1.5 pH

The slight pH increase is caused mainly by a larger availability of inorganic carbon and thus

hydrogencabonae anion

Figure 26: pH

3.1.6 Gas Phase Response

The gas signal increases with the sludge increase, while the methane fraction decreases due

to a more pronounced carbon dioxide generation. The extremely low hydrogen content,

Results and Discussion

73

which should be at around 1% showcases another problem with the model. Either the mass

transfer coefficient must be increased or the uptake by the hydrogenotrophes must be

decreased. Figure 27 shows the trend of gas flowrate increase, while Figure 28 underlines

the actual methane ratio diminishment due to a consistent CO2 increase.

Figure 27: gas flowrate signal in response to sludge modulation

Figure 28: Gas fractions in the outflow as response to sludge modulation signal

3.1.7 Mass Balances

The execution of mass balances of COD, nitrogen and carbon showed positive results with

regards to COD (Figure 29), while showcasing problems with nitrogen and carbon. It is

evident from Figure 30 and Figure 31 that in both cases the output does not converge to the

Results and Discussion

74

input Such discrepancies are due to bad parametric characterizations and should be

investigated further to improve the model correctness.

With regard to the COD balance, the Figure 32 helps the reader to better understand

the COD mass flows.

Figure 29: COD mass balance

Figure 30: Nitrogen Mass Balance

Results and Discussion

75

Figure 31: Carbon Mass Balance

Figure 32: Sankey Diagram for the COD flow

3.2 PARAMETER SENSITIVITY: VFA ACCUMULATION

The extreme soluble substratum consumption prompted a biochemical parameters

sensitivity testing with regards to the acetoclastic methanogens. The goal was to simulate

pH drop due to the VFA accumulation and to obtain an output of VFA concentration

comparable to that of the dataset. To obtain this, it was decided to act on the biomass’

dynamic by changing the decay rate and the half-saturation constant. The parameter

tweaking was performed with inhibition dynamic disabled at first, then with inhibition

enabled. The inhibition-free dynamic garnered positive results, obtaining the target

concentration of 1 [g/l] for the parameter set expressed in Table 20.

Results and Discussion

76

KS,X2 KD,2 S2

[mgCOD/l] [1/d] [mgCOD/l]

0.04 0.02 5.7

0.15 0.02 16.32

0.15 0.2 29.2

0.15 0.25 1000

Table 20: Parameter Tweaking Results

The resulting VFA output is closer to the real data, as shown in the Figure 33

Figure 33: Comparison between the real VFA output(blue) and the model output(blue) [mgCOD/l]

Results and Discussion

77

The increased VFA output influences the sludge pH, dropping it form the default steady

state value of 6.9 to 6.6:

Figure 34: pH signal with no inhibition

The repetition of the simulation with the inhibition dynamic on the other hand, provided

negative results: there was no parameter combination that could provide a sufficiently high

VFA output without crashing the system. The acetoclasts’ growth rate becomes reasonably

inhibited by pH values below 6.7 (Rosen, 2006) so an accumulation of VFA leads to an

positive feedback. The pH drop leads to a rapid inorganic carbon stripping, showcased in

Figure 37, which crashes the system because the inorganic carbon concentration becomes

close to zero, and inorganic carbon induced inhibition diverging to infinity. Since the uptake

of VFA is strongly inhibited, its accumulation occurs, as shown in Figure 35

Figure 35: VFA output with the inhibition dynamic present.

Results and Discussion

78

The pH dynamic has physical significance only until the day three, when the non-negative

conditions are violated.

Figure 36: pH with inhibition present

Figure 37: Inorganic Carbon stripping

Results and Discussion

79

3.3 HYDROGEN INJECTION IMPLEMENTATION

A preliminary hydrogen feed response was performed, which aimed at the evaluation of

the system response to high levels of gaseous hydrogen. The hydrogen injection was

implemented by regulating a Fick mass transfer equation with a parameter regulated by a PI

controller. At the pre-hydrogen bubbling phase, the gaseous carbon dioxide flowrate was

obtained:

pre-hydro CO2 flowrate

qco2 1369.35 [ml/d]

ṅco2 0.054 [mmol/d]

This parameter was then used to obtain the hydrogen flowrates at different stages:

H2/CO2 0 0.76 1.01 1.9 2.6 4.36 6.01 7

starting day

0 21 41 50 59 68 95 108 [d]

qh2 0 1040.71 1383.04 2601.77 3560.31 5970.37 8229.79 9585.45 [ml/d]

ṅh2 0 56.03 74.46 140.06 191.67 321.41 443.05 516.03 [mmol/d]

Table 21: Hydrogen bubbling Values

The target signal and the response obtained by the controller until the simulation crash are

showcased in the figure below.

Figure 38: Hydrogen Inflow

Results and Discussion

80

The results of the simulation once again underlined the main problem of the model, that the

substrate consumption is too rapid. The hydrogenotrophic methanogens use the dissolved

inorganic carbon as a carbon source, so an increase in hydrogen availability leads to a higher

carbon uptake.

3.3.1 On exceeding the stoichiometric ratio

The hydrogen stress test evidenced a carbon depletion dynamic. For hydrogen

injection ratios higher than the stoichiometric ratio (which is 4), the carbon uptake becomes

too excessive with respect to the supply. This leads to a greater degree of inorganic carbon

consumption as observed in Figure 39 where a pronounced decrease of inorganic carbon

availability (which is comprised mainly by hydrogencarbonate) can be observed in

concomitance with the beginning of the fifth hydrogen injection phase at around the day 70.

The consumption of all available inorganic carbon eventually reaches values close zero.

Figure 39: Inorganic Carbon

The lack of inorganic carbon inhibits the growth of hydrogenotrophic fauna (X3),

event that leads to the accumulation of hydrogen. The Hydrogen on the other hand inhibits

the acetogens and acidogens (X1). The decay of X1 means that, except for the input flow,

there is no substrate available for the acetoclasts. This eventually leads for the biomass

washout, as shown in the figure below. It is observable that the pronounced biomass

Results and Discussion

81

diminishment occurs when the inorganic carbon is exhausted, at around the day 90 (Figure

40).

Figure 40: Biomass Washout

The Upgrading of the biogas peaked at 99 [%] when the H2/CO2 ratio reached the

stoichiometric value of 4. To mimic more closely the data values, a new parameter was

introduced, which function is to partition the influent hydrogen between the aqueous and

gaseous phase. The comparison between the scenarios without and with the bypass

parameter can be seen in Figure 41 and Figure 42.

Figure 41: Output Gas Fractioning - no bypass

Results and Discussion

82

Figure 42: Influent Gas Fractioning with bypass parameter

This has generated more acceptable results, as shown in the table below. The CO2 is still

very low, but this is mainly due to high carbon uptake on X3’s part.

Max

Dataset

Value

fbypass = 0 fbypass = 0.05

CO2 16.7 1 1

CH4 79.2 99 92

H2 3.9 0 7

[%] [%] [%]

Table 22: Gas Fractioning -: Comparison between the data and the model

Results and Discussion

83

3.4 PH CONTROL IMPLEMENTATION

Rather that modeling UIT’s data for the acid/base dosage, the pH control was implemented

by acting upon the cation and anion compartments by the means of a PI controller, whose

target signal was the average dataset’s pH of 7.3. The pH control starts at day 20,

concurrently with the beginning of hydrogen bubbling.

Figure 43: pH signal comparison

Results and Discussion

84

3.5 DATA FIT

To judge how well the model simulates the real data, the comparison between the

model’s output and the data was made. Although the experimentation time lasts for 140

days, the simulation stops at approximately 120 days because of numerical error – namely

because inhibition functions (equations 2.32 - 2.33) diverge to infinity as inorganic carbon

reaches values close to zero. In this paragraph the model’s output is confronted with the

real data.

3.5.1 Total Solids

The model’s total solid output closely resembles that of the data. The fluctuations

occur because the feed sludge is stocked in tanks that have different concentrations and

that are changed throughout the experiment.

Figure 44: Total Solids in the outflow - data vs model

Figure 45: Total Solids in the inflow – data vs model

Results and Discussion

85

3.5.2 Soluble COD

The soluble COD comparison once again evidenced the presence of a substrate uptake rate

that is too high – the model’s output is three orders of magnitude lower than that of the

data. The model’s data skyrockets at day 100 because the substrate uptake is arrested by

the inhibitory dynamics, therefore leading to the substrate accumulation.

Figure 46: Total Soluble COD

3.5.3 Alkalinity

The alkalinity mimics pretty closely the data until the hydrogen injection reaches levels

beyond the stoichiometric (day 68): the discrepancy between the carbon consumption of

model’s and data’s hydrogenotrophic methanogens is once again evident.

Figure 47: Alkalinity - Comparison between the model and data

Results and Discussion

86

3.5.4 Gas Phase

The model’s gas output mimics pretty closely that of the data, diverging for higher degrees

of hydrogen injection and for when the initial conditions still have a relevant impact. This

shows the validity of the Fick transfer parameters and gas dynamic alike.

Figure 48: Total Gas Outflow - model vs data

Figure 49: methane and carbon dioxide gas outflow

The gas fractioning on the other hand, shows a greater degree of upgrading for the model’s

data, mainly because the model’s biomass has a higher degree of substrate uptake than the

Results and Discussion

87

real biomass, leading on one hand to a greater degree of methane production and on the

other hand to a higher inorganic carbon consumption.

Figure 50: Methane and Carbon Dioxide fractioning in the output gas

The model’s hydrogen gas outflow resembles less closely the data but at least, after the

introduction of the bypass parameter it is not null.

Figure 51: Hydrogen Gas Outflow Comparison

Results and Discussion

88

3.6 LFT PARAMETER IDENTIFICATION

The main goal of the LFT implementation is the identification of key parameters within

the AMOCO model. To perform a feasibility of such endeavor, a preemptive phase of

validation of the LFT model was performed. At first the data from the Modelica simulation

was obtained, which was then implemented in the MATLAB algorithm. The LFT algorithm

was then tweaked by changing

• sets of output variables,

• simulation times,

• number of iterations

• interpolation methods

• LFT linearization algorithms of the DAE system

Since the LFT model, as proposed by (Rosito, 2018) didn’t seem to work, a wast array of

models with different ways of linearizing were developed. The most successful LFT form as

well as parametric setup that delivered the most precise results are reported in the code

section. The results that stemmed from the conclusion of the validation phase led to the

following conclusions.

• The precision depended upon the simulation time

• The simulation time longer that 10000 [s] sensibly dropping down the precision.

Further increments in the simulation time crashed the algorithm execution.

• The parametric identification was the most precise when the simulation time was

kept at around 1000 [s]

Results and Discussion

89

The results of the comparison between the real and obtained parameters for each

simulation time are shown in the table below. The mean squared error (MSE) is reported at

the bottom of the table.

PARAMETER TRUE

VALUES

100

[S]

1000

[S]

10000

[S]

100000

[S]

200000

[S]

FVFA 0.48225 0.99487 0.40909 0.39354 1.00000 ERROR

KLA,H2 0.00049 0.05739 0.00036 0.00029 2.51350 ERROR

KLA 0.00231 0.00249 0.00654 0.00140 0.04109 ERROR

MSE 0.00000 0.08867 0.00179 0.00262 2.19493 N/A

Table 23: Error as a function of time

The following graph showcases the error dynamic as a function of the simulation time.

Figure 52: Error Dynamic as a function of simulation time

0

0.5

1

1.5

2

2.5

10 100 1000 10000 100000

[s]

Error Dynamic

Mean Squared Error

Chapter 4

4 CONCLUSIONS

The aim of this work is the validation of the modified AMOCO model and the feasibility of

the identification of key parameters using LFT-based algorithm.

The original model developed by (Rosito, 2019) was deemed too incomplete and

therefore was subjected to an extensive review. A better ionic dynamic was added and the

pH equation was made useable. The ionic and hydrogen control of the original model were

formulated without an actual way to implement them, so a way was developed by using PI

controllers. More than that, a review of the parameters was made, underlining and

correcting incongruencies. The LFT model as formulated originally by (Rosito, 2019) was

binned because did not deliver acceptable results, and a new model was developed.

With regards to the model’s validation it was concluded that a more precise

parametric identification of the substrate uptake is necessary. Since the consumption of

soluble substrate was deemed too rapid an attempt to identify a parametric set for the

growth dynamic of acetoclastic bacteria was made. Although such attempt was fruitful in

the absence of inhibition, in its presence no such parametric combination was obtained. It is

therefore advised to identify a new parametric set for the acidogens & acetogens (X1) and

acetoclasts (X2) that works well with the hydrogen and inorganic carbon-induced inhibition

constants.

Despite the initial premises, the Fick mass transfer coefficients seemed to perform

quite well to forecast the gas dynamic for low levels of hydrogen bubbling. The main issue

with gaseous phase dynamic is the hydrogen uptake rate. The hydrogen injection stress test

evidenced that the hydrogenotrophic activity consumes all the available hydrogen. This

behavior leads to gaseous hydrogen levels that are not in line with the data. To overcome

this without modifying the uptake parameters, a new parameter was introduced. Rather

than inferring about a new Fick mass transfer constant KLA,H2, for the bubbled hydrogen, it

Conclusions

91

was decided to model the injected H2 bypass to the reactor’s headspace. This harnessed

acceptable yet questionable results: the possibility of implementation of both bypass

parameter and Fick constant for the bubbled hydrogen should be explored.

Although the COD mass balance provided positive results, the carbon and nitrogen

mass balances could not achieve enough convergence. A review of carbon and nitrogen

partitioning parameters therefore is advised, starting with volatile solids’ partitioning

parameters Cx01, Cnd, Nx01, Nnd.

On a positive note it must be said that the model behaves reasonably well for the

upgrading with the H2/CO2 ratio below the stoichiometric. The problems arise when the

carbon induced inhibition starts to take effect, when the inorganic carbon levels drop too

low.

With regards to the implementation of the LFT algorithm for the parametric identification,

mildly successful results were obtained. Although the simulation time could not be

extended beyond one day, an acceptable result was obtained for a simulation time of 1000s.

Ways of simplifying the LFT formulation should be explored. Possible results could be

obtained by simulating the model in days instead of seconds.

BIBLIOGRAPHY

Panpan Li, Chao He, Gang Li, Pan Ding, Mingming Lan, Zan Gao, and Youzhou Jiao (2020).

Biological pretreatment of corn straw for enhancing degradation efficiency and biogas

production.

Scotto di Perta, E., Cervelli, E., Pironti di Campagna, M., Pindozzi, S (2020). From biogas to

biomethane: Techno-economic analysis of an anaerobic digestion power plant in a

cattle/buffalo farm in central Italy.

E.Ryckebosch, M.Drouillon, H.Vervaeren (2011). Techniques for transformation of biogas to

biomethane.

Chen, Xiao & Vinh, Hoang & Avalos Ramirez, Antonio & Rodrigue, Denis & Kaliaguine, Serge.

(2015). Membrane gas separation technologies for biogas upgrading.

J. Patrick Murphy, Arthur Wellinger, David Baxter (2013). The Biogas Handbook: Science,

Production and Applications

J. D. Murphy, T. Thamisisroj (2003). Fundamental science and engineering of the anaerobic

digestion process for biogas production

Bibliography

93

Buswell AM and Hatfield WD (1936) Bulletin No. 32, Anaerobic Fermentations. State of

Illinois Department of Registration and Education.

SUEZ degremont® water handbook (2007).

C. Rosén & U. Jeppsson (2006). Aspects on ADM1 Implementation within the BSM2

Framework.

S. Civardi (2013). Identificazione Parametrica di Modelli di Digestione Anaerobica

Xiao Yuan Chen, Hoang Vinh-Thang, Antonio Avalos Ramirez, Denis Rodriguez (2015).

Membrane gas separation technologies for biogas upgrading.

F. Malpei, D. Gardoni (2008). La digestione anaerobica e i principi del processo biologico e di

dimensionamento

M. Rosito (2019). Biological biogas upgrading by hydrogenotrophic methanogenesis:

optimization and modeling.

A. Santus (2019). Ex situ biological biogas upgrading model implementation.

A. Negri (2017). Modellizzazione e controllo di un reattore a biogas.

R. Speece (1997). Anaerobic biotechnology for industrial wastewater

G. Tchobanglous (2006). Wastewater engineering: treatment and reuse.

Bibliography

94

A. Wellinger, J. Murphy (2013). The biogas handbook: science production and application.

Decreto interministeriale 2 marzo 2018 - Promozione dell’uso del biometano nel settore dei

trasporti

G. Ferretti A. Dalla Bona F. Malpei E. Ficara (2015). LFT modelling and identification of

anaerobic digestion.

U. Bünger et al., 2014. Power-to-Gas (PtG) in transport Status quo and perspectives for

development

BP Statical Review of World Energy 67th Edition [Report]. - 2018.

APPENDIX A

A.1 VARIABLES

VARIABLE DESCRIPTION MEASURE

X01 Degradable Particulate Matter kgCOD/m3

Xnd Non-degradable Particulate

Matter

kgCOD/m3

X1 Acidogenic Flora kgCOD/m3

X2 Acetoclasic Methanogens kgCOD/m3

X3 Hydrogenotrophic Methanogens kgCOD/m3

S1 Soluble COD kgCOD/m3

S2 VFAS kgCOD/m3

S4 Dissolved Hydrogen kgCOD/m3

S5 Dissolved Methane kgCOD/m3

IN Inorganic Nitrogen M

Snh3 Ammonia M

IC Inorganic Carbon M

Shco3 Hydrogencarbonate M

Scat Cations M

San Anions M

Sgas,h2 Gaseous Hydrogen M

Sgas,co2 Gaseous Carbon Dioxide M

Sgas,ch4 Gaseous Methane M

S3 Dissolved CO2 M

Snh4 Ammonium M

Appendix A

96

A.2 PETERSEN MATRIX

rate

𝐾ℎ1𝑋01

𝜇1 max

𝑌 x1

𝑠 1𝑆 1+𝐾𝑠1

I 𝑝𝐻1I 𝐻2X1

Y𝑥2

𝜇2 max

𝑌x2

𝑠 2

𝑆2+𝐾𝑠2

I 𝑝𝐻2X2

Y𝑥3

𝜇3 max

𝑌x3

𝑠 4

𝑆4+𝐾𝑠42

I 𝑝𝐻3I 𝐼𝐶X3

𝐾d1𝑋1

𝐾d2𝑋2

𝐾d3𝑋3

𝑘𝐿𝑎𝐻2( 𝑆4−16𝐾𝐻,𝐻2𝑃 𝑔𝑎𝑠,𝐻2)

𝑘𝐿𝑎( 𝑆5/64−𝐾𝐻,𝐶𝐻5𝑃 𝑔𝑎𝑠,𝐶𝐻4)

𝑘𝐿𝑎( 𝑆3−𝐾𝐻,𝐶𝑂2𝑃 𝑔𝑎𝑠,𝐶𝑂2)

𝑘𝐿𝑎𝐻2( 𝐾𝐻,𝐻2𝑃𝑏𝑢𝑏𝑏𝑙𝑒,𝐻2

−𝑆 4/64) (−VALVh2)

IN

Nx0

1 -

Ns1

Ns1

- Y

x1N

bac

-Yx2

Nb

ac

-Yx3

Nb

ac

Nb

ac -

Nx0

1

Nb

ac -

Nx0

1

Nb

ac -

Nx0

1

IC

Cx0

1-C

s1

(Cs1

- (Y

x1C

bac

+

(1-

Yx1

)Cs2

fvfa

))

(Cs2

- (

Y x2C

bac

+ (1

- Y

x2)

Cch

4))

(-(Y

x3C

bac

+ (1

-

Y x3)

Cch

4))

Cb

ac-C

x01

Cb

ac-C

x01

Cb

ac-C

x01

-1

S 5

1 -

Yx2

1 -

Yx3

-64

S 4

(1 -

Yx1

)fh

2

-1

-16

+16

(1-

f byp

ass)

S 2

(1 -

Yx1

)

f vfa

-1

S 1

1-k

nd

-1

X3

Y x3

-1

X2

Y x2

-1

X1

Y x1

-1

Xn

d

k nd

X0

1

-1

1

1

1

Pro

cess

Dis

inte

grat

io

n a

nd

Hyd

roly

sis

S 1 u

pta

ke

S 2 u

pta

ke

S 4 u

pta

ke

X1

dec

ay

X2

dec

ay

X3

dec

ay

H2

gas

tran

sfer

CH

4 m

ass

tran

sfer

CO

2 m

ass

tran

sfer

Hyd

roge

n

inje

ctio

n

Appendix A

97

A.3 MODEL’S EQUATIONS

𝑑𝑋01𝑑𝑡

=𝑄

𝑉(𝑋01,𝑖𝑛 − 𝑋01) + 𝑘ℎ1𝑋01 + 𝑘𝑑,1𝑋1 + 𝑘𝑑,2𝑋2 + 𝑘𝑑,3𝑋3 A.3.1

𝑑𝑋𝑛𝑑𝑑𝑡

=𝑄

𝑉(−𝑋𝑛𝑑) + 𝑘𝑛𝑑𝑘ℎ1𝑋01 A.3.2

𝑑𝑋1𝑑𝑡

=𝑄

𝑉(𝑋1,𝑖𝑛 − 𝑋1) + 𝑌1

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2 − 𝑘𝑑,1𝑋1 A.3.3

𝑑𝑋2𝑑𝑡

=𝑄

𝑉(𝑋2,𝑖𝑛 − 𝑋2) + 𝑌2

𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2 − 𝑘𝑑,2𝑋2 A.3.4

𝑑𝑋3𝑑𝑡

=𝑄

𝑉(𝑋3,𝑖𝑛 − 𝑋3) + 𝑌3

𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶 − 𝑘𝑑,3𝑋3 A.3.5

𝑑𝑆1𝑑𝑡=𝑄

𝑉(𝑆1,𝑖𝑛 − 𝑆1) + (1 − 𝑘𝑛𝑑)𝑘ℎ1𝑋01 −

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2 A.3.6

𝑑𝑆2𝑑𝑡=𝑄

𝑉(𝑆2,𝑖𝑛 − 𝑆2) + (1 − 𝑌𝑥,1)𝑓𝑉𝐹𝐴

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2

−𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2

A.3.7

𝑑𝑆4𝑑𝑡=𝑄

𝑉(𝑆4,𝑖𝑛 − 𝑆4) + (1 − 𝑌𝑥,1)𝑓𝐻2

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2

−𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶−16𝑘𝑙𝑎,𝐻2 (𝑆416− 𝑘𝐻,ℎ2𝑝ℎ2)

+(1 − 𝑓𝑏𝑦𝑝𝑎𝑠𝑠)16𝑘𝑙𝑎,𝐻2(𝑆416− 𝑘𝐻,ℎ2𝑝𝑏𝑢𝑏𝑏𝑙𝑒)𝑉𝐴𝐿𝑉𝐻2

A.3.8

𝑑𝑆5𝑑𝑡=𝑄

𝑉(𝑆5,𝑖𝑛 − 𝑆5) + (1 − 𝑌𝑥,1)

𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2

+(1 − 𝑌𝑥,3)𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶

−64𝑘𝑙𝑎 (𝑆564− 𝑘𝐻,𝑐ℎ4𝑝𝑐ℎ4)

A.3.9

𝑑𝐼𝐶

𝑑𝑡=𝑄

𝑉(𝐼𝐶𝑖𝑛 − 𝐼𝐶) + (𝐶𝑆1 − (1 − 𝑌𝑥,1)𝑓𝑉𝐹𝐴 − 𝑌𝑥,1𝐶𝑏𝑎𝑐)

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2

+(𝐶𝑠2 − (1 − 𝑌𝑥,2)𝐶𝑐ℎ4 − 𝑌𝑥,2𝐶𝑏𝑎𝑐)𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2

+(−(1 − 𝑌𝑥,3)𝐶𝑐ℎ4 − 𝑌𝑥,3𝐶𝑏𝑎𝑐)𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶−𝑘𝑙𝑎(𝑆3 − 𝑘𝐻,𝑐𝑜2𝑝𝑐𝑜2)

+(𝐶𝑥01 − 𝐶𝑥1)𝑘ℎ1𝑋01

A.3.10

Appendix A

98

𝑑𝑆ℎ𝑐𝑜3𝑑𝑡

= 𝐾𝑎𝑏,𝑐𝑜2(𝑆ℎ𝑐𝑜3 − (𝑘𝑎,𝑐𝑜2 + 𝑆ℎ) − 𝑘𝑎,𝑐𝑜2𝐼𝐶) A.3.11

𝐼𝐶 = 𝑆3 + 𝑆ℎ𝑐𝑜3 A.3.12

𝑑𝐼𝑁

𝑑𝑡=𝑄

𝑉(𝐼𝑁𝑖𝑛 − 𝐼𝑁) + (𝑁𝑥01 − 𝑁𝑥1)𝑘ℎ1𝑋01 + (𝑁𝑠1 − 𝑌𝑥1𝑁𝑏𝑎𝑐)

𝜇𝑚𝑎𝑥,1𝑌𝑥,1

𝑆1𝐾𝑠,1 + 𝑆1

𝑋1𝐼𝑝𝐻,1𝐼𝐻2

−𝑌𝑥,2𝑁𝑏𝑎𝑐𝜇𝑚𝑎𝑥,2𝑌𝑥,2

𝑆2𝐾𝑠,2 + 𝑆2

𝑋2𝐼𝑝𝐻,2

−𝑌𝑥,3𝑁𝑏𝑎𝑐𝜇𝑚𝑎𝑥,3𝑌𝑥,3

𝑆4𝐾𝑠,4 + 𝑆4

𝑋3𝐼𝑝𝐻,3𝐼𝐼𝐶 + (𝑁𝑏𝑎𝑐 + 𝑁𝑥01)(𝑘𝑑,1𝑋1 + 𝑘𝑑,2𝑋2 + 𝑘𝑑,3𝑋3)

A.3.13

𝑑𝑆𝑛ℎ3𝑑𝑡

= 𝐾𝑎𝑏,𝐼𝑁(𝑆𝑛ℎ3 − (𝑘𝑎,𝐼𝑁 + 𝑆ℎ) − 𝑘𝑎,𝐼𝑁𝐼𝑁)

A.3.14

𝐼𝑁 = 𝑆𝑛ℎ3 + 𝑆𝑛ℎ4 A.3.15

𝑑𝑆𝑐𝑎𝑡𝑑𝑡

=𝑄

𝑉(𝑆𝑐𝑎𝑡,𝑖𝑛 − 𝑆𝑐𝑎𝑡) +

𝑞𝑎𝐶𝑎𝑉𝑉𝐴𝐿𝑉𝑎 A.3.16

𝑑𝑆𝑎𝑛𝑑𝑡

=𝑄

𝑉(𝑆𝑎𝑛,𝑖𝑛 − 𝑆𝑎𝑛) +

𝑞𝑏𝐶𝑏𝑉𝑉𝐴𝐿𝑉𝑏 A.3.17

𝑑𝑆𝑔𝑎𝑠,ℎ2

𝑑𝑡= −𝑞𝑔𝑎𝑠

𝑆𝑔𝑎𝑠,ℎ2

𝑉𝑔𝑎𝑠+𝑉

𝑉𝑔𝑎𝑠𝑘𝑙𝑎,ℎ2(

𝑆416− 𝑘𝐻,ℎ2𝑝ℎ2)

+ 𝑓𝑏𝑦𝑝𝑎𝑠𝑠𝑘𝑙𝑎,ℎ2(𝑘𝐻,ℎ2𝑝ℎ2,𝑏𝑢𝑏𝑏𝑙𝑒 −𝑆416)

A.3.18

𝑑𝑆𝑔𝑎𝑠,𝑐ℎ4

𝑑𝑡= −𝑞𝑔𝑎𝑠

𝑆𝑔𝑎𝑠,𝑐ℎ4

𝑉𝑔𝑎𝑠+𝑉

𝑉𝑔𝑎𝑠𝑘𝑙𝑎 (

𝑆564− 𝑘𝐻,𝑐ℎ4𝑝𝑐ℎ4) A.3.19

𝑑𝑆𝑔𝑎𝑠,𝑐𝑜2

𝑑𝑡= −𝑞𝑔𝑎𝑠

𝑆𝑔𝑎𝑠,𝑐𝑜2

𝑉𝑔𝑎𝑠+𝑉

𝑉𝑔𝑎𝑠𝑘𝑙𝑎(𝑆3 − 𝑘𝐻,𝑐𝑜2𝑝𝑐𝑜2)

A.3.20

𝑞𝑔𝑎𝑠 =𝑅𝑇

(𝑝𝑔𝑎𝑠 − 𝑝𝑣)𝑉(𝑘𝑙𝑎,ℎ2(

𝑆416− 𝑘𝐻,ℎ2𝑝ℎ2) + 𝑘𝑙𝑎 (

𝑆564− 𝑘𝐻,𝑐ℎ4𝑝𝑐ℎ4) + 𝑘𝑙𝑎(𝑆3

− 𝑘𝐻,𝑐𝑜2𝑝𝑐𝑜2))

A.3.21

𝜃 = 𝑆𝑐𝑎𝑡 + 𝑆𝑛ℎ4 − 𝑆ℎ𝑐𝑜3 −𝑆216− 𝑆𝑎𝑛 A.3.22

𝑆ℎ =−𝜃 + √𝜃2 − 4𝑘𝑤

2

A.3.23

𝑝𝐻 = −𝑙𝑜𝑔(𝑆ℎ) A.3.24

Appendix A

99

IpH1 = {exp(−3(

pH − pHUL,1pHUL,1 − pHLL,1

)

2

) ∶ pH < pHUL.1

1 ∶ pH > pHUL,1

A.3.25

IpH2 = {exp(−3(

pH − pHUL,2rHUL,2 − pHLL,2

)

2

) ∶ pH < pHUL,2

1 ∶ pH > pHUL,2

A.3.26

lpH3 = {exp(−3(

pH − pHUL,3pHUL,3 − pHLL.3

)

2

) ∶ pH < pHUL,3

1 ∶ pH > pHUL,3

A.3.27

IH2 =1

1 +S4KI,H2

A.3.28

𝐼𝐼𝐶 =1

1 +𝐾𝐼,𝐼𝐶𝐼𝑐

A.3.29

APPENDIX B

OPENMODELICA CODE

The code is structured in a series of nested packages. The first four packages contain records

of declaration of types, parameters, variables and data. The main is comprised by

AMOCO_vanilla and its ancillary ControlRoom model. AMOCO_vanilla contains the record

invocation and the AMOCO equations and does not work standalone. The ControlRoom

provides a way to control the model by modulating the COD input, pH, hydrogen injection.

The target signals for pH and hydrogen injection are also constructed here by invoking the

two blocks – Picontroller and StepSignal. DataSim deals with the creation of time-dependent

variables that describe the dataset’s inputs and outputs. DataConversion converts the

model’s variables from a strictly SI unit format to a more readable one. The COD mass

balance is also calculated here. Finally, CarbonNitrogenBalance is used to perform the

carbon and nitrogen mass balance.

The model is extensively commented to guide the programmer in its use and

understanding.

Note on OpenModelica version

It is advisable to utilize the OpenModelica version 12/19 or earlier since the latest version

changes the connector syntax and does not work with the presented model.

Appendix B

101

package AMOCO

package Types

// this package contains the declaration of types

type Conc_COD = Real(unit = "kgCOD/m3", min = 0);

type Conc_COD_model = Real(unit = "mgCOD/l", min = 0);

type M = Real(unit = "kmol/m3", min = 0);

type M_model = Real(unit = "mmol/l", min = 0) "molarity in

mg/l";

type Rate_COD = Real(unit = "kgCOD/(m3*s)");

type Rate_COD_Mass = Real(unit = "kgCOD/s", min = 0);

type Rate_M = Real(unit = "kmol/(m3*s)");

type Rate = Real(unit = "1/s");

type Pressure = Real(unit = "Pa", min = 0);

type Temp = Real(unit = "K", min = 0);

type Flow = Real(unit = "m3/s", min = 0);

type Volume = Real(unit = "m3", min = 0);

type Ratio_COD = Real(unit = "kgCOD/kgCOD", min = 0);

type Henry = Real(unit = "kmol/(m3*Pa)", min = 0);

type MassFlow = Real(unit = "kg/s", min = 0);

type MassFlowDay = Real(unit = "g/d", min = 0);

end Types;

package Parameters

// this package contains the declarations of model's

parameters.

// The contestual èarameters are grouped in records

// Each record is invoked by other objects

record BioChemicalParameters

//biochemical parameters

parameter Real R(unit = "kJ/(kmol*K)") =

Modelica.Constants.R;

parameter AMOCO.Types.Rate mu_1_max = 3.5 / 86400

"changed from 3.5";

parameter AMOCO.Types.Rate mu_2_max = 0.33 / 86400

"default = 0.33";

parameter AMOCO.Types.Rate mu_3_max = 2 / 86400 "default

= 2";

parameter AMOCO.Types.Rate k_d_1 = 0.02 / 86400 "default

= 0.02";

parameter AMOCO.Types.Rate k_d_2 = 0.02 / 86400 "default

= 0.02, ";

parameter AMOCO.Types.Rate k_d_3 = 0.02 / 86400 "default

= 0.02";

parameter AMOCO.Types.Conc_COD K_s_1 = 0.5 "default =

0.5";

parameter AMOCO.Types.Conc_COD K_s_2 = 0.04 "changed

from 0.04, threshold = 0.151725";

parameter AMOCO.Types.Conc_COD K_s_4 = 7e-6 "default =

7e-6";

parameter AMOCO.Types.Ratio_COD Y_x_1 = 0.15 "default =

0.15";

Appendix B

102

parameter AMOCO.Types.Ratio_COD Y_x_2 = 0.05 "default =

0.05";

parameter AMOCO.Types.Ratio_COD Y_x_3 = 0.06 "default =

0.06";

parameter AMOCO.Types.Rate k_h1 = 0.3855 / 86400 "";

parameter AMOCO.Types.Ratio_COD k_nd = 0.599 "";

end BioChemicalParameters;

record PhysioChemicalParameters

// Physiochemical parameters. This record invokes the

// Reactor parameters and biochemical parameters for

// it's values computation

extends AMOCO.Parameters.ReactorParameters;

extends AMOCO.Parameters.BioChemicalParameters;

parameter AMOCO.Types.Rate k_la = 8.3 / 3600;

parameter AMOCO.Types.Rate k_la_h2 = 1.75 / 3600;

parameter AMOCO.Types.Pressure p_h2o = 0.0313 * exp(5290

* (1 / 298 - 1 / T)) * 1e5;

parameter AMOCO.Types.Henry k_h_co2 = 0.035 * exp(-19410

/ R * (1 / 298 - 1 / T)) / 1e5;

parameter AMOCO.Types.Henry k_h_ch4 = 0.0014 * exp(-

14240 / R * (1 / 298 - 1 / T)) / 1e5;

parameter AMOCO.Types.Henry k_h_h2 = 0.00078 * exp(-4180

/ R * (1 / 298 - 1 / T)) / 1e5;

parameter Real k_w = exp(148.9802 - 13847.26 / T -

23.6521 * log(T));

// mol^2/L^2 Kw

parameter Real k_ab_co2(unit = "1/(M*s)") = 115740;

//rosen:1e10/86400

parameter Real k_ab_IN(unit = "1/(M*s)") = 115740;

parameter AMOCO.Types.M k_a_co2 = 10 ^ (-6.35) *

exp(7646 / R * (1 / 298 - 1 / T));

parameter AMOCO.Types.M k_a_IN = 10 ^ (-9.25) *

exp(51965 / R * (1 / 298 - 1 / T));

parameter AMOCO.Types.M k_a_ac = 10 ^ (-4.86);

end PhysioChemicalParameters;

record ReactorParameters

// parameters describing the reactor conditions

AMOCO.Types.Flow Q = 0.0005 / 86400;

parameter AMOCO.Types.Volume V = 11e-3;

//parameter Real HRT(unit="d") = V/Q/86400;

parameter AMOCO.Types.Volume V_gas = 3e-3;

parameter AMOCO.Types.Temp T = 310.15;

parameter AMOCO.Types.Pressure p_gas = 101325;

end ReactorParameters;

record PartiotioningParameters " nitrogen and carbon

partitioning parameters "

parameter Real C_bac(unit = "kmolC/kgCOD") = 0.0313;

parameter Real C_x01(unit = "kmolC/kgCOD") = 0.02786;

Appendix B

103

parameter Real C_s1(unit = "kmolC/kgCOD") = 0.0313;

parameter Real C_s2(unit = "kmolC/kgCOD") = 0.0313;

parameter Real C_ch4(unit = "kmolC/kgCOD") = 0.0156;

parameter Real N_bac(unit = "kmolN/kgCOD") = 0.00625;

parameter Real N_x01(unit = "kmolN/kgCOD") = 0.002;

parameter Real N_s1(unit = "kmolN/kgCOD") = 0.007 * 0.3;

parameter Real N_aa(unit = "kmolN/kgCOD") = 0.007 * 0.3

"aminoacid's Nitrogen fraction";

parameter Real f_bypass = 0.05 "frazione di idrogeno

gorgogliato trasferito direttamente in fase gassosa";

end PartiotioningParameters;

record CODFractioningParameters

// parameters describing how the COD is subdivided between

// hydrogen and VFA

parameter Real f_h2 = f_h2_acido + f_h2_aceto;

parameter Real f_h2_acido = 0.19 * 0.47 + 0.06 * 0.265 +

0.3 * 0.265;

parameter Real f_h2_aceto = (0.13 * 0.20 + 0.43 * 0.27 +

0.15 * 0.23) * 0.81 + 0.19;

parameter Real f_vfa = 1 - f_h2;

end CODFractioningParameters;

record InputParameters

// Input parameters for an old version of the model

// based on the works of ROSito

extends AMOCO.Parameters.PhysioChemicalParameters;

parameter AMOCO.Types.Conc_COD X_01_in = 18;

parameter AMOCO.Types.Conc_COD X_1_in = 0.01;

parameter AMOCO.Types.Conc_COD X_2_in = 0.01;

parameter AMOCO.Types.Conc_COD X_3_in = 0.01;

parameter AMOCO.Types.Conc_COD S_1_in = 4.7 - 0.48078;

parameter AMOCO.Types.Conc_COD S_2_in = 0.48078;

parameter AMOCO.Types.Conc_COD S_4_in = 0;

parameter AMOCO.Types.Conc_COD S_5_in = 1e-5;

parameter Real Alk_in(unit = "mgCaCO3/m3") = 1500;

parameter AMOCO.Types.M IN_in = 0.014061308;

parameter AMOCO.Types.M S_nh4_in = 124 / 1020 * 1e-3;

parameter Real pH_in = 6.67;

parameter AMOCO.Types.M S_cat_in = 0.02;

parameter AMOCO.Types.M S_an_in = 0.011355;

parameter AMOCO.Types.M Alk_in_hco3eq = Alk_in / 50 /

1000;

parameter AMOCO.Types.M S_2_ion_in = k_a_ac * S_2_in /

64 / (k_a_ac + 10 ^ (-pH_in));

parameter AMOCO.Types.M OH_in = k_w / 10 ^ (-pH_in);

parameter AMOCO.Types.M IC_in = Alk_in_hco3eq - OH_in -

S_2_ion_in - S_nh4_in + 10 ^ (-pH_in);

parameter AMOCO.Types.M HCO3_in = k_a_co2 * IC_in /

(k_a_co2 + 10 ^ (-pH_in));

end InputParameters;

Appendix B

104

record InputParameters2020

// input paramters based on the dataset's values

extends AMOCO.Parameters.PhysioChemicalParameters;

//from reactor 2 data (phase 0)

parameter AMOCO.Types.Conc_COD X_01_in = 26.75690669

"Dataset averaged on pre upgrade phase";

parameter AMOCO.Types.Conc_COD X_1_in = 0.01 "rosen";

parameter AMOCO.Types.Conc_COD X_2_in = 0.01 "rosen";

parameter AMOCO.Types.Conc_COD X_3_in = 0.01 "rosen";

parameter AMOCO.Types.Conc_COD S_1_in = (4.7 -

0.48078)/18*X_01_in "pesato sui valori di bresso assay";

parameter AMOCO.Types.Conc_COD S_2_in =

0.48078/18*X_01_in "pesato sui valori di bresso assay";

parameter AMOCO.Types.Conc_COD S_4_in = 0;

parameter AMOCO.Types.Conc_COD S_5_in = 0;

parameter Real Alk_in(unit = "gCaCO3/m3") = 1600

"dataset";

parameter AMOCO.Types.M IN_in = S_nh4_in*(10^(-

pH_in)+k_a_IN)/10^(-pH_in) "";

parameter AMOCO.Types.M S_nh4_in = 253/18.04/1000

"bresso assay 2";

parameter Real pH_in = 6.67;

parameter AMOCO.Types.M S_cat_in = 0.02 "Rosen";

parameter AMOCO.Types.M S_an_in = 10^(-pH_in) + S_nh4_in

+ S_cat_in - S_2_ion_in - k_w/(10^(-pH_in)) - HCO3_in

"0.011355";

parameter AMOCO.Types.M Alk_in_hco3eq = Alk_in / 50 /

1000;

parameter AMOCO.Types.M S_2_ion_in = k_a_ac * S_2_in /

64 / (k_a_ac + 10 ^ (-pH_in));

parameter AMOCO.Types.M OH_in = k_w / 10 ^ (-pH_in);

parameter AMOCO.Types.M IC_in = Alk_in_hco3eq - OH_in -

S_2_ion_in - S_nh4_in + 10 ^ (-pH_in);

parameter AMOCO.Types.M HCO3_in = k_a_co2 * IC_in /

(k_a_co2 + 10 ^ (-pH_in));

end InputParameters2020;

record InhibitionParameters

parameter Real pH_ll_1 = 4.5;

parameter Real pH_ul_1 = 5.5;

parameter Real pH_ll_2 = 5.8;

parameter Real pH_ul_2 = 6.7;

parameter Real pH_ll_3 = 5;

parameter Real pH_ul_3 = 6;

parameter AMOCO.Types.M K_I_IC = 1e-4;

parameter AMOCO.Types.Conc_COD K_I_h2 = 1e-5;

end InhibitionParameters;

record InitialConditions

//initial concentration values of the state variables

Appendix B

105

// as calculated by rosito

parameter AMOCO.Types.Conc_COD X_01_start = 18;

parameter AMOCO.Types.Conc_COD X_nd_start = 1.1e-2;

parameter AMOCO.Types.Conc_COD X_1_start = 0.3;

parameter AMOCO.Types.Conc_COD X_2_start = 0.6;

parameter AMOCO.Types.Conc_COD X_3_start = 0.4;

parameter AMOCO.Types.Conc_COD S_1_start = 4.7 -

0.48078;

parameter AMOCO.Types.Conc_COD S_2_start = 0.48078;

parameter AMOCO.Types.Conc_COD S_4_start = 1e-5;

parameter AMOCO.Types.Conc_COD S_5_start = 1e-5;

parameter AMOCO.Types.M IN_start = 0.014061308;

parameter AMOCO.Types.M S_nh4_start = 0.01402439;

parameter AMOCO.Types.M S_nh3_start = 3.69183e-5;

parameter AMOCO.Types.M IC_start = 0.022565578;

parameter AMOCO.Types.M S_hco3_start = 0.01526711;

parameter AMOCO.Types.M S_3_start = 0.007298468;

parameter AMOCO.Types.M S_cat_start = 0.02;

parameter AMOCO.Types.M S_an_start = 0.0114;

parameter AMOCO.Types.M S_oh_start = 1.1215e-7;

parameter AMOCO.Types.M S_gas_h2_start = 1e-7;

parameter AMOCO.Types.M S_gas_ch4_start = 1e-7;

parameter AMOCO.Types.M S_gas_co2_start = 1e-7;

parameter Real pH_start = 6.67;

end InitialConditions;

record InitialConditions2020

// initial concentration values of the state variables

// as obtained drom the dataset's data

extends AMOCO.Parameters.PhysioChemicalParameters;

parameter AMOCO.Types.Conc_COD X_01_start = 16.79786517-

(X_1_start + X_2_start + X_3_start + X_nd_start) "Bresso";

parameter AMOCO.Types.Conc_COD X_nd_start = 1.1e-2;

parameter AMOCO.Types.Conc_COD X_1_start = 0.3 "Rosen";

parameter AMOCO.Types.Conc_COD X_2_start = 0.6 "Rosen";

parameter AMOCO.Types.Conc_COD X_3_start = 0.4 "Rosen";

parameter AMOCO.Types.Conc_COD S_1_start = (2234.277778

- 1070.434579)/1000 "Bresso";

parameter AMOCO.Types.Conc_COD S_2_start =

1070.434579/1000 "Bresso";

parameter AMOCO.Types.Conc_COD S_4_start = 1e-5 "Rosen";

parameter AMOCO.Types.Conc_COD S_5_start = 1e-5 "Rosen";

parameter AMOCO.Types.M IN_start = S_nh4_start*(10^(-

pH_start)+k_a_IN)/10^(-pH_start);

parameter AMOCO.Types.M S_nh4_start = 253/18.04/1000

"Bresso";

parameter AMOCO.Types.M S_nh3_start = IN_start-

S_nh4_start;

//parameter AMOCO.Types.M IC_start = 0.058169616;

parameter AMOCO.Types.M S_hco3_start = 0.05312533;

parameter AMOCO.Types.M S_3_start = 0.0050443;

Appendix B

106

parameter AMOCO.Types.M S_oh_start = k_w/10^(-pH_start);

parameter AMOCO.Types.M S_cat_start = 0.02;

parameter AMOCO.Types.M S_an_start = 10^(-pH_start) +

S_nh4_start + S_cat_start - S_2_ion_start - k_w/(10^(-

pH_start)) - S_hco3_start "0.0114";

parameter AMOCO.Types.M S_gas_h2_start = 1e-7;

parameter AMOCO.Types.M S_gas_ch4_start = 1e-7;

parameter AMOCO.Types.M S_gas_co2_start = 1e-7;

parameter Real pH_start = 7.32;

parameter Real Alk_start(unit = "mgCaCO3/l") = 4446;

parameter AMOCO.Types.M Alk_eq_start = Alk_start / 50 /

1000;

parameter AMOCO.Types.M S_2_ion_start = k_a_ac *

S_2_start / 64 / (k_a_ac + 10 ^ (-pH_start));

parameter AMOCO.Types.M IC_start = Alk_eq_start -

S_oh_start - S_2_ion_start - S_nh4_start + 10 ^ (-pH_start)

"0.058169616";

end InitialConditions2020;

end Parameters;

package Variables

// this package contains the declaration of variables

record StateVariables

// defition of the state variables and the invocation

// of initial conditions

// si puo invocare due parametri iniziali:

//InitialConditions - rosito

//InitialConditions2020 - condizioni dataset

extends AMOCO.Parameters.InitialConditions2020;

AMOCO.Types.Conc_COD X_01(start = X_01_start, fixed =

true);

AMOCO.Types.Conc_COD X_nd(start = X_nd_start, fixed =

true);

AMOCO.Types.Conc_COD X_1(start = X_1_start, fixed =

true);

AMOCO.Types.Conc_COD X_2(start = X_2_start, fixed =

true);

AMOCO.Types.Conc_COD X_3(start = X_3_start, fixed =

true);

AMOCO.Types.Conc_COD S_1(start = S_1_start, fixed =

true);

AMOCO.Types.Conc_COD S_2(start = S_2_start, fixed =

true);

AMOCO.Types.Conc_COD S_4(start = S_4_start, fixed =

true);

AMOCO.Types.Conc_COD S_5(start = S_5_start, fixed =

true);

AMOCO.Types.M S_gas_co2(start = S_gas_co2_start, fixed =

true);

Appendix B

107

AMOCO.Types.M S_gas_h2(start = S_gas_h2_start, fixed =

true);

AMOCO.Types.M S_gas_ch4(start = S_gas_ch4_start, fixed =

true);

AMOCO.Types.M IN(start = IN_start, fixed = true);

AMOCO.Types.M S_nh3(start = S_nh3_start, fixed = false);

AMOCO.Types.M S_nh4(start = S_nh4_start, fixed = true);

AMOCO.Types.M IC(start = IC_start, fixed = true);

AMOCO.Types.M S_hco3(start = S_hco3_start, fixed =

true);

AMOCO.Types.M S_3(start = S_3_start, fixed = false);

AMOCO.Types.M S_an(start = S_an_start, fixed = false);

AMOCO.Types.M S_cat(start = S_cat_start, fixed = true);

AMOCO.Types.M S_h(start = 10 ^ (-pH_start), fixed =

true);

Real pH(start = pH_start, fixed = false);

end StateVariables;

record AlgebraicVariables

//AMOCO.Types.M S_nh4(start = 0.01402439, fixed =

false);

//AMOCO.Types.M S_h(start = 10 ^ (-6.67), fixed =

false);

//AMOCO.Types.M S_3(start = 0.007298468, fixed = false);

/* bilancio di cariche */

AMOCO.Types.M S_2_ion;

Real Alk(unit = "mgCaCO3/l");

Real THETA;

end AlgebraicVariables;

record RateVariables " variables used in rate equations "

AMOCO.Types.Rate_COD r_1;

AMOCO.Types.Rate_COD r_2;

AMOCO.Types.Rate_COD r_3;

AMOCO.Types.Rate_COD r_4;

AMOCO.Types.Rate_COD r_5;

AMOCO.Types.Rate_COD r_6;

AMOCO.Types.Rate_COD r_7;

AMOCO.Types.Rate_M r_8;

AMOCO.Types.Rate_M r_9;

AMOCO.Types.Rate_M r_10;

AMOCO.Types.Rate_M r_11;

AMOCO.Types.Rate_M r_12;

AMOCO.Types.Rate_M r_13;

end RateVariables;

record InhibitionVariables

Real I_pH_1;

Real I_pH_2;

Real I_pH_3;

Real I_IC;

Appendix B

108

Real I_h2;

end InhibitionVariables;

record GasVariables

AMOCO.Types.Pressure p_h2;

AMOCO.Types.Pressure p_ch4;

AMOCO.Types.Pressure p_co2;

AMOCO.Types.Pressure p_tot;

AMOCO.Types.Flow q_gas;

Real x_h2;

Real x_ch4;

Real x_co2;

end GasVariables;

end Variables;

package Dataset2020

// this package contains the time series of the data

obtained from

// the dataset. The missing data was interpolated by spline

method

record TimeSeries

// definition of the time vecttor spanning through 140

days

parameter Real TimeSeries[:, 1] = [0; 86400; 172800;

259200; 345600; 432000; 518400; 604800; 691200; 777600;

864000; 950400; 1036800; 1123200; 1209600; 1296000; 1382400;

1468800; 1555200; 1641600; 1728000; 1814400; 1900800; 1987200;

2073600; 2160000; 2246400; 2332800; 2419200; 2505600; 2592000;

2678400; 2764800; 2851200; 2937600; 3024000; 3110400; 3196800;

3283200; 3369600; 3456000; 3542400; 3628800; 3715200; 3801600;

3888000; 3974400; 4060800; 4147200; 4233600; 4320000; 4406400;

4492800; 4579200; 4665600; 4752000; 4838400; 4924800; 5011200;

5097600; 5184000; 5270400; 5356800; 5443200; 5529600; 5616000;

5702400; 5788800; 5875200; 5961600; 6048000; 6134400; 6220800;

6307200; 6393600; 6480000; 6566400; 6652800; 6739200; 6825600;

6912000; 6998400; 7084800; 7171200; 7257600; 7344000; 7430400;

7516800; 7603200; 7689600; 7776000; 7862400; 7948800; 8035200;

8121600; 8208000; 8294400; 8380800; 8467200; 8553600; 8640000;

8726400; 8812800; 8899200; 8985600; 9072000; 9158400; 9244800;

9331200; 9417600; 9504000; 9590400; 9676800; 9763200; 9849600;

9936000; 10022400; 10108800; 10195200; 10281600; 10368000;

10454400; 10540800; 10627200; 10713600; 10800000; 10886400;

10972800; 11059200; 11145600; 11232000; 11318400; 11404800;

11491200; 11577600; 11664000; 11750400; 11836800; 11923200;

12009600; 12096000; 12182400];

end TimeSeries;

record Alk

//definition of matrix containing output alkalinity and

time vector

//invocation of record containiing the time series

Appendix B

109

extends AMOCO.Dataset2020.TimeSeries ;

parameter Real DataFrameAlk[:, 2] = [TimeSeries,

DataAlk];

//parameter Real DataSeriesAMOCO[:,1]

=fill(1,size(TimeSeriesAMOCO,1),1);

parameter Real DataAlk[:, 1] = [4446; 4490.527873; 4467;

4406.222127; 4339; 4289.862352; 4258.229863; 4237.246197;

4220.055021; 4199.8; 4180.040664; 4206; 4308.702839;

4420.71793; 4450; 4341.980778; 4192; 4106.632337; 4087.48418;

4109.919856; 4149.303688; 4181; 4182.024696; 4136;

4041.337509; 3949; 3910.01608; 3923.126489; 3974; 4028; 4026;

4067.780948; 4119; 4117.684942; 4082.201526; 4051;

4051.375159; 4066; 4074.008454; 4085; 4109.189324;

4128.786395; 4119; 4068.274932; 4018; 4007.387542; 4014;

4007.023303; 3983.787238; 3948.657554; 3906; 3860.5834;

3818.788872; 3787.400609; 3773.202803; 3782.979647;

3823.515333; 3901.594053; 4024; 4195.463573; 4412.5;

4660.079313; 4885.205942; 5025.392914; 5018.153257; 4801;

4367.69557; 3937; 3724.450839; 3707.698114; 3804.91997;

3934.29455; 4014; 3991.199953; 3929; 3900.993705; 3906.78457;

3927.478583; 3944.181731; 3938; 3898.488647; 3849; 3815.79661;

3802.984832; 3809.131974; 3832.805344; 3872.57225; 3927;

3993.498025; 4064.844248; 4132.658718; 4188.56148;

4224.172584; 4231.112074; 4201; 4130.762115; 4038.546999;

3947.808942; 3882.002231; 3864.581154; 3919; 4044.94123; 4147;

4134.480653; 4030.611578; 3887.102176; 3755.66185; 3688;

3715.094861; 3785; 3838.054105; 3866.661022; 3876.240887;

3872.213834; 3860; 3843.682953; 3822; 3792.454901;

3752.963494; 3701.544636; 3636.217186; 3555; 3465.040455;

3410; 3429.370739; 3509.453321; 3623.250535; 3743.765165;

3844; 3899.752647; 3898; 3834.276436; 3727.166973;

3601.019293; 3480.181075; 3389; 3351.823748; 3393;

3536.876436; 3807.800735; 4230.120579; 4828.183648];

end Alk;

record VFA

// output vfa data

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrameVFA[:, 2] = [TimeSeries,

DataVFA];

parameter Real DataVFA[:, 1] = [1070.434579;

841.1805447; 678.7226055; 577.1167289; 530.4188819;

532.6850317; 577.9711454; 660.3331899; 773.8271323;

912.5089398; 1070.434579; 1241.660018; 1420.241223;

1600.234162; 1775.6948; 1940.679107; 2089.243047; 2215.44259;

2313.333701; 2376.972348; 2400.414497; 2377.716117;

2302.933173; 2170.121634; 1973.337466; 1706.636636;

1385.138559; 1108.216444; 996.3069489; 1169.846729;

1681.760836; 2314.927903; 2784.714953; 2873.997152;

2635.682214; 2190.185988; 1657.924328; 1159.313084;

800.5474749; 630.9401869; 663.5873524; 823.6014186;

Appendix B

110

1014.098911; 1138.196357; 1099.01028; 852.8719955;

568.9719626; 423.2413573; 405.7150669; 459.9539065;

529.5186916; 572.3566252; 603.9604623; 654.2093458;

742.2300715; 844.1400469; 925.3043323; 951.0879882;

886.8560748; 723.8794719; 557.0523364; 480.7552887;

483.6915239; 528.144881; 576.3991993; 590.7383178;

544.4236971; 454.6272849; 349.4986503; 257.1873624;

205.8429905; 223.6151036; 338.653271; 549.1923576;

733.8084112; 775.5604387; 695.0981006; 547.4687207;

387.7196229; 270.8981308; 239.7341167; 287.6876454; 395.90133;

545.5177837; 717.6796195; 893.5294507; 1054.20989;

1180.863551; 1258.772141; 1289.773739; 1279.845522;

1234.964663; 1161.108337; 1064.253719; 950.3779824;

825.4583031; 695.4718551; 566.3958132; 444.2073518;

334.8836456; 244.4018692; 179.1459844; 147.1271028;

153.2776487; 188.5881492; 240.5636572; 296.7092254;

344.5299065; 373.4990951; 380.9635514; 368.1156227;

343.6566368; 318.1651666; 302.2197852; 306.3990654;

336.562992; 379.6971963; 420.1002373; 450.196739; 464.4428418;

457.2946858; 423.2084112; 362.7768064; 301.1392523;

263.3382041; 249.481821; 253.4446886; 269.1013922;

290.3265173; 310.9946494; 324.9803738; 327.7175418;

320.8770675; 307.6891306; 291.3839111; 275.1915888;

262.3423435; 256.0663551; 259.5938035; 276.1548684;

308.9797298; 361.2985675];

end VFA;

record CodSol

//output soluble cod data

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrameCodSol[:, 2] = [TimeSeries,

DataCodSol];

parameter Real DataCodSol[:, 1] = [2234.277778;

2203.732832; 2234.277778; 2312.371539; 2424.473042;

2557.041212; 2696.534974; 2829.413252; 2942.134974;

3021.159063; 3052.944444; 3028.817383; 2959.573499;

2860.875749; 2748.387092; 2637.770486; 2544.688889;

2482.706757; 2456.994533; 2470.62416; 2526.667579;

2628.196731; 2778.283557; 2980; 3207.591813; 3320;

3203.087928; 2965.713162; 2772.481815; 2788; 3107.351132;

3547.527827; 3856; 3847.659646; 3607.087078; 3286.284686;

3037.254863; 3012; 3296.919399; 3716; 4066.703609;

4302.803842; 4417.152271; 4402.600467; 4252; 3977.445021;

3668; 3414.167916; 3235.235469; 3132.685287; 3108;

3156.293355; 3247.203575; 3344; 3415.850355; 3455.515895;

3461.656257; 3432.931079; 3368; 3268.254066; 3146.009958;

3016.315768; 2894.219586; 2794.769504; 2733.013612; 2724;

2778.141348; 2887.308689; 3038.737643; 3219.66383;

3417.322873; 3618.95039; 3811.782003; 3983.053333; 4120;

4212.072744; 4257.582782; 4257.056447; 4211.020075; 4120;

3990.383298; 3852; 3734.095116; 3640.130599; 3567.122639;

Appendix B

111

3512.087426; 3472.04115; 3444; 3424.743387; 3410.103602;

3395.676157; 3377.056562; 3349.84033; 3309.622972; 3252;

3175.114043; 3087.296203; 2999.424702; 2922.37776;

2867.033598; 2844.270437; 2864.966497; 2940; 3071.393431;

3225.746335; 3360.802524; 3434.305809; 3404; 3256.011865;

3090; 3006.957644; 2999.684788; 3033.933111; 3075.454288;

3090; 3059.425456; 3030; 3051.031229; 3119.565633;

3219.584422; 3335.068807; 3450; 3550.343255; 3630;

3685.803943; 3718.381743; 3729.308297; 3720.158503;

3692.507257; 3647.929457; 3588; 3515.148189; 3435.22095;

3354.919617; 3280.945523; 3220; 3178.784381; 3164;

3182.348189; 3240.53028; 3345.247608; 3503.201504];

end CodSol;

record Q_gas

//gas outflow

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrameQ_gas[:, 2] = [TimeSeries,

DataQ_gas];

parameter Real DataQ_gas[:, 1] = [

3497.503513 ;

3352.529274 ;

3406.894614 ;

3524.686183 ;

3497.503513 ;

3524.686183 ;

3488.442623 ;

3606.234192 ;

3705.903981 ;

3397.833724 ;

3252.859485 ;

4729.784543 ;

4222.374707 ;

4104.583138 ;

3923.36534 ;

3669.660422 ;

3678.721311 ;

3515.625293 ;

3678.721311 ;

3859.93911 ;

3742.147541 ;

4045.468798 ;

3696.281374 ;

4040.267351 ;

4029.919358 ;

4004.234939 ;

3882.47704 ;

4049.217281 ;

4044.207055 ;

4068.583536 ;

4184.978717 ;

Appendix B

112

4546.257585 ;

4487.259692 ;

4186.095513 ;

4364.212234 ;

4520.303196 ;

4381.884248 ;

4240.08628 ;

4041.380957 ;

3972.757454 ;

4113.057011 ;

4231.071555 ;

4350.416828 ;

4150.677753 ;

4222.208397 ;

4258.212078 ;

4269.419906 ;

4271.491437 ;

4280.446346 ;

4357.085706 ;

4308.837774 ;

4059.574608 ;

4482.48598 ;

4583.109784 ;

4628.112471 ;

4627.289809 ;

4751.579368 ;

4710.318055 ;

4800.610592 ;

5115.06178 ;

5333.85315 ;

5382.446014 ;

5293.830199 ;

5223.288204 ;

5101.865562 ;

5366.461667 ;

6014.400901 ;

5792.928538 ;

5514.011886 ;

5413.562406 ;

5403.907803 ;

5397.375782 ;

5306.294045 ;

5093.129451 ;

4920.905467 ;

4950.47296 ;

5133.435792 ;

5369.086074 ;

5556.715915 ;

5595.617427 ;

5441.162226 ;

5273.039956 ;

Appendix B

113

5277.386614 ;

5441.805568 ;

5704.267031 ;

6002.741215 ;

6275.198333 ;

6459.608597 ;

6509.386627 ;

6439.724671 ;

6281.259385 ;

6064.627423 ;

5820.465442 ;

5579.410096 ;

5372.098042 ;

5222.287895 ;

5126.226116 ;

5073.281126 ;

5052.821346 ;

5054.215197 ;

5066.831102 ;

5090.96514 ;

5170.624032 ;

5342.774756 ;

5572.514678 ;

5806.973759 ;

5993.281962 ;

6078.569248 ;

6044.337865 ;

6009.579196 ;

6076.782334 ;

6144.938075 ;

6062.162639 ;

5676.572248 ;

4921.215959 ;

4068.874176 ;

3477.260138 ;

3410.3456 ;

3757.136386 ;

4312.896831 ;

4872.891275 ;

5232.384055 ;

5261.155046 ;

5127.046274 ;

5030.511743 ;

5004.391225 ;

5039.620934 ;

5127.137086 ;

5257.875893 ;

5422.773571 ;

5612.766333 ;

5812.720825 ;

5983.225411 ;

Appendix B

114

6078.798889 ;

6053.960054 ;

5863.227702 ;

5461.12063 ;

4802.157634 ;

4802.157634 ;

4802.157634 ;

4802.157634 ;

4802.157634

];

end Q_gas;

record x_co2

// co2 fraction

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrame_x_co2[:, 2] = [TimeSeries,

Data_x_co2];

parameter Real Data_x_co2[:, 1] = [23; 22.99102921;

22.7; 22.47552156; 22.3689847; 22.49659249; 22.4; 22.22967338;

22.14367177; 22.2; 21.868637; 21.8494788; 22.37921; 23.08108;

23.65771905; 23.92307158; 24; 23.90945758; 23.95215577;

23.95820584; 23.92557331; 24; 24; 23.88075052; 23.85312935;

23.94652295; 23.96251739; 23.84617524; 23.68332175;

23.59054242; 23.4; 23.20751043; 23.10180807; 23.1;

23.27260083; 23.2171647; 22.86926287; 22.87454798;

23.01905424; 22.74489229; 22.26418637; 22.14089013; 22.2;

22.23273106; 22.24979138; 22.08045897; 21.98866481; 21.9;

21.9; 21.9; 21.71070932; 21.61341209; 21.45090403;

21.32162726; 21.3; 21.21276287; 21.19930459; 21.08191933;

20.86817234; 20.6; 20.8; 20.8; 20.9; 20.8; 20.6; 20.7; 20.7;

20.7; 20.93661947; 21.26186883; 21.42880847; 21.19049874;

20.3; 18.688538; 17; 15.90345452; 15.45603371; 15.56188565;

16.12515839; 17.05; 18.19043578; 19.2; 19.76636285;

19.9142292; 19.75256256; 19.39032646; 18.93648443; 18.5;

18.16860183; 17.94507913; 17.81098625; 17.74787755;

17.73730737; 17.76083007; 17.8; 17.83484069; 17.83925246;

17.78560478; 17.64626715; 17.39360906; 17; 16.4528836; 15.8;

15.10538291; 14.43794363; 13.86781291; 13.46512146; 13.3;

13.40970552; 13.7; 14.04807117; 14.34830476; 14.49938596;

14.4; 13.97898425; 13.28578483; 12.4; 11.41908217;

10.51190023; 9.865177214; 9.665636129; 10.1; 11.19032501;

12.3; 12.80059581; 12.75507866; 12.39926361; 11.9689657; 11.7;

11.74201588; 11.9; 11.94806792; 11.88151347; 11.75092506;

11.6068911; 11.5; 11.48084016; 11.6; 11.6; 11.6; 11.6; 11.6];

end x_co2;

record x_ch4

// methane fraction

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrame_x_ch4[:, 2] = [TimeSeries,

Data_x_ch4];

Appendix B

115

parameter Real Data_x_ch4[:, 1] = [71.28554552; 71.6;

71.6; 71.43184979; 71.3689847; 71.5; 71.5; 71.5; 71.71835883;

72; 72.1988178; 72.25017373; 71.82079; 71.26182; 70.84228095;

70.55399583; 70.4; 69.76620306; 69.43922114; 69.28358832;

69.47671994; 69.7105007; 69.8; 69.68075052; 69.22809458;

68.76043115; 68.63748261; 68.7; 68.75003475; 69.00945758;

69.2; 69.68087622; 69.89638387; 69.99770515; 70.17260083;

70.2828353; 70.38268428; 70.62364395; 70.87858136;

70.79708131; 70.78108484; 70.95910987; 70.97579972; 71.1;

71.07489569; 71.03908206; 71.18866481; 71.1; 71.1;

71.00041725; 70.2433936; 70.54635163; 70.6; 70.44325452;

70.31912378; 70.03828861; 70.09095967; 70.10897079;

70.21591383; 70.3; 70.1; 69.9; 69.8; 69.9; 69.8; 69.9; 70.2;

69.6; 68.77864794; 68.05313494; 67.49829797; 67.188974; 67.2;

67.53340672; 67.9; 68.03049344; 68.0264566; 68.08217305;

68.39192633; 69.15; 70.40214728; 71.6; 72.1956383;

72.23705685; 71.92122898; 71.445128; 71.00572723; 70.8;

70.97463237; 71.4751614; 72.19683691; 73.03490871;

73.88462662; 74.64124044; 75.2; 75.48934211; 75.57045159;

75.53770029; 75.48546003; 75.50810266; 75.7; 76.10939312;

76.6; 77.02631351; 77.39044913; 77.73142799; 78.08827124;

78.5; 78.95011059; 79.2; 79.04905933; 78.67075534;

78.33207367; 78.3; 78.71880578; 79.24190564; 79.4;

78.85827006; 77.81982006; 76.62223504; 75.60310001; 75.1;

75.27648939; 75.6; 75.51359773; 75.05900723; 74.42761787;

73.81081901; 73.4; 73.35381783; 73.7; 74.40121198;

75.29080069; 76.16978341; 76.83917742; 77.1; 76.75326843;

75.6; 75.6; 75.6; 75.6; 75.6];

end x_ch4;

record pH

// ph output

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFramepH[:, 2] = [TimeSeries, DatapH];

parameter Real DatapH[:, 1] = [7.32; 7.32; 7.32; 7.31;

7.31; 7.31; 7.31; 7.31; 7.32; 7.32; 7.32; 7.31; 7.28; 7.23;

7.22; 7.23; 7.23; 7.23; 7.21; 7.21; 7.22; 7.22; 7.22; 7.23;

7.24; 7.22; 7.22; 7.21; 7.22; 7.22; 7.22; 7.23; 7.22; 7.22;

7.22; 7.22; 7.22; 7.23; 7.23; 7.22; 7.23; 7.23; 7.22; 7.22;

7.22; 7.22; 7.23; 7.22; 7.22; 7.22; 7.22; 7.23; 7.22; 7.21;

7.22; 7.22; 7.22; 7.23; 7.23; 7.23; 7.22; 7.23; 7.31; 7.23;

7.23; 7.23; 7.23; 7.23; 7.28; 7.33; 7.33; 7.35; 7.37; 7.37;

7.36; 7.35; 7.32; 7.30; 7.31; 7.33; 7.33; 7.33; 7.34; 7.35;

7.35; 7.38; 7.45; 7.41; 7.39; 7.38; 7.37; 7.36; 7.36; 7.36;

7.37; 7.38; 7.39; 7.39; 7.38; 7.39; 7.42; 7.44; 7.44; 7.45;

7.49; 7.54; 7.47; 7.45; 7.42; 7.40; 7.40; 7.39; 7.42; 7.49;

7.54; 7.52; 7.73; 7.83; 7.89; 7.64; 7.36; 7.37; 7.42; 7.45;

7.46; 7.41; 7.35; 7.32; 7.36; 7.39; 7.41; 7.43; 7.46; 7.56;

7.56; 7.58; 7.39; 7.44; 7.50; 7.55; 7.45; 7.40];

end pH;

Appendix B

116

record x_h2

//hydrogen fraction

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrame_x_h2[:, 2] = [TimeSeries,

Data_x_h2];

parameter Real Data_x_h2[:, 1] = [0.0001; 0.43326418;

0.815766831; 1.150072844; 1.438647115; 1.683954535;

1.888459999; 2.054628399; 2.184924628; 2.28181358;

2.347760149; 2.385229226; 2.396685706; 2.384594482;

2.351420447; 2.299628494; 2.231683517; 2.150050408;

2.057194061; 1.95557937; 1.847671227; 1.735934525;

1.622834158; 1.51083502; 1.402402003; 1.3; 1.206093905;

1.123148611; 1.053629012; 1; 0.962537379; 0.932760594; 0.9;

0.855944905; 0.801720432; 0.740810654; 0.676699647;

0.612871486; 0.552810245; 0.5; 0.457223832; 0.424460849;

0.400989165; 0.386086893; 0.379032149; 0.379103046;

0.385577698; 0.39773422; 0.414850726; 0.436205329;

0.461076145; 0.488741287; 0.518478869; 0.549567005;

0.58128381; 0.612907398; 0.643715883; 0.672987379; 0.7;

0.718486726; 0.7; 0.628441367; 0.55928932; 0.56591659;

0.721695907; 1.1; 1.744632703; 2.581122261; 3.505428019;

4.413509325; 5.201325526; 5.764835969; 6; 5.86913619; 5.6;

5.435712137; 5.415417392; 5.527266578; 5.75941051; 6.1;

6.504795998; 6.8; 6.826538251; 6.61379545; 6.238270912;

5.776463949; 5.304873874; 4.9; 4.620382811; 4.452727467;

4.365780299; 4.328287638; 4.308995814; 4.276651158; 4.2;

4.060476398; 3.890265316; 3.734239445; 3.637271475;

3.644234096; 3.8; 4.107262825; 4.4; 4.515139948; 4.470131273;

4.327552624; 4.14998265; 4; 3.923544415; 3.9; 3.893833234;

3.876395722; 3.820760349; 3.7; 3.539865415; 3.576818746; 4.1;

5.27751538; 6.793335879; 8.210398689; 9.091640999; 9;

7.767185779; 6.3; 5.532468366; 5.432421143; 5.726139736;

6.139905553; 6.4; 6.299666932; 5.9; 5.322706336; 4.664097206;

4.014134908; 3.46278174; 3.1; 3.015751987; 3.3; 3.3; 3.3; 3.3;

3.3];

end x_h2;

record COD_x_out

//particulate COD outflow

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrame_COD_x_out[:, 2] = [TimeSeries,

Data_COD_x_out];

parameter Real Data_COD_x_out[:, 1] = [16.79786517;

16.82178904; 16.79786517; 16.74230994; 16.67133968;

16.60117077; 16.54801954; 16.52810236; 16.55763558;

16.65283555; 16.82991864; 17.09004265; 17.3741313;

17.60804975; 17.71766317; 17.62883674; 17.26743562;

16.60748591; 15.81565742; 15.10678091; 14.69568711;

14.79720678; 15.62617065; 17.2664457; 19.27804387;

21.09001331; 22.13140219; 21.83125866; 19.6186309;

16.06262333; 14.98191697; 16.09265503; 17.80033225;

Appendix B

117

18.81780665; 19.08738945; 18.85875521; 18.38157844;

17.9055337; 17.62385998; 17.50405417; 17.4687961; 17.4872395;

17.54015656; 17.60831948; 17.67250046; 17.76046439;

18.08794701; 18.8488546; 19.96180496; 21.27659378;

22.64301673; 23.91086948; 24.92994772; 25.55004712;

25.62096336; 24.9924921; 23.51442904; 21.03656985; 17.4087102;

12.88083444; 9.30368166; 8.395966625; 9.747552224;

12.41608837; 15.45922497; 17.93461194; 19.1754682;

19.61728868; 19.85140632; 19.99023005; 20.0264378;

19.95270748; 19.76171702; 19.47267687; 19.21092762;

19.09565777; 19.11531748; 19.22567233; 19.38248788;

19.54152971; 19.66679321; 19.7551931; 19.81187394;

19.84198028; 19.85065665; 19.84304762; 19.82429772;

19.79955152; 19.77298238; 19.74487901; 19.71455894;

19.68133968; 19.64453877; 19.60347375; 19.55746212;

19.50888641; 19.47238903; 19.46567737; 19.50645884;

19.61244084; 19.80133075; 20.05300701; 20.1960321;

20.07655911; 19.7624193; 19.37686352; 19.04314261;

18.88450739; 18.95152891; 19.00405895; 18.81035726;

18.4630347; 18.13578986; 18.00232136; 18.2363278; 18.89791308;

19.5928022; 19.91664962; 19.8792063; 19.59374735; 19.17354789;

18.73188304; 18.39331431; 18.32754883; 18.65385212;

19.24457749; 19.91035022; 20.46179557; 20.70953882;

20.55309401; 20.24753014; 20.0681799; 20.01587546;

20.02282391; 20.02123235; 19.94330783; 19.72125746;

19.2872883; 19.2872883; 19.2872883; 19.2872883; 19.2872883];

end COD_x_out;

record COD_x_in

//particulate COD inflow

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrame_COD_x_in[:, 2] = [TimeSeries,

Data_COD_x_in];

parameter Real Data_COD_x_in[:, 1] = [24.38097614;

24.38097614; 24.38097614; 24.38097614; 24.38097614;

26.92090577; 26.92090577; 26.92090577; 26.92090577;

29.46083541; 29.46083541; 29.46083541; 27.81382794;

27.81382794; 27.81382794; 26.16682048; 26.16682048;

26.16682048; 27.32736176; 27.32736176; 27.32736176;

27.32736176; 27.32736176; 27.32736176; 28.92334912;

28.92334912; 28.92334912; 28.92334912; 28.92334912;

29.24127183; 29.24127183; 29.24127183; 29.52592639;

30.94919918; 30.94919918; 30.94919918; 30.72190238;

30.72190238; 30.72190238; 29.91056299; 29.91056299;

29.91056299; 29.91056299; 29.67849576; 29.67849576;

29.76583283; 29.76583283; 30.19520522; 30.19520522;

30.19520522; 30.19520522; 30.19520522; 30.19520522;

32.59960341; 32.59960341; 32.59960341; 32.59960341;

32.59960341; 32.59960341; 32.59960341; 32.59960341;

32.59960341; 32.59960341; 32.59960341; 32.59960341;

32.59960341; 32.59960341; 32.59960341; 29.94671618;

Appendix B

118

29.94671618; 29.94671618; 29.94671618; 29.94671618;

29.94671618; 29.94671618; 29.94671618; 29.94671618;

29.94671618; 29.94671618; 27.29382896; 27.29382896;

27.29382896; 27.29382896; 27.29382896; 27.29382896;

27.29382896; 27.29382896; 28.50072332; 28.50072332;

28.50072332; 28.50072332; 28.50072332; 28.50072332;

28.50072332; 29.70761769; 29.70761769; 29.70761769;

29.70761769; 28.07092842; 28.07092842; 28.07092842;

28.07092842; 28.07092842; 36.35200495; 36.35200495;

36.35200495; 36.35200495; 36.35200495; 36.35200495;

18.46763937; 18.46763937; 18.46763937; 18.46763937;

18.46763937; 18.46763937; 18.46763937; 18.46763937;

18.46763937; 18.46763937; 30.01043718; 30.01043718;

30.01043718; 30.01043718; 30.01043718; 30.01043718;

30.01043718; 30.01043718; 30.01043718; 29.75734841;

29.75734841; 29.75734841; 29.75734841; 29.75734841;

29.75734841; 29.75734841; 14.17144359; 14.17144359;

14.17144359; 14.17144359; 14.17144359; 14.17144359; 14.17144];

end COD_x_in;

record VFA_in

//vfa inflow

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrame_VFA_in[:, 2] = [TimeSeries,

Data_VFA_in];

parameter Real Data_VFA_in[:, 1] = [0.651215873;

0.651215873; 0.651215873; 0.651215873; 0.651215873;

0.719057393; 0.719057393; 0.719057393; 0.719057393;

0.786898914; 0.786898914; 0.786898914; 0.742907344;

0.742907344; 0.742907344; 0.698915775; 0.698915775;

0.698915775; 0.729913833; 0.729913833; 0.729913833;

0.729913833; 0.729913833; 0.729913833; 0.772542655;

0.772542655; 0.772542655; 0.772542655; 0.772542655;

0.781034371; 0.781034371; 0.781034371; 0.788637494;

0.82665311; 0.82665311; 0.82665311; 0.820582013; 0.820582013;

0.820582013; 0.798911137; 0.798911137; 0.798911137;

0.798911137; 0.792712622; 0.792712622; 0.795045395;

0.795045395; 0.806513932; 0.806513932; 0.806513932;

0.806513932; 0.806513932; 0.806513932; 0.870735407;

0.870735407; 0.870735407; 0.870735407; 0.870735407;

0.870735407; 0.870735407; 0.870735407; 0.870735407;

0.870735407; 0.870735407; 0.870735407; 0.870735407;

0.870735407; 0.870735407; 0.799876789; 0.799876789;

0.799876789; 0.799876789; 0.799876789; 0.799876789;

0.799876789; 0.799876789; 0.799876789; 0.799876789;

0.799876789; 0.729018172; 0.729018172; 0.729018172;

0.729018172; 0.729018172; 0.729018172; 0.729018172;

0.729018172; 0.76125432; 0.76125432; 0.76125432; 0.76125432;

0.76125432; 0.76125432; 0.76125432; 0.793490468; 0.793490468;

0.793490468; 0.793490468; 0.749774498; 0.749774498;

0.749774498; 0.749774498; 0.749774498; 0.970962052;

Appendix B

119

0.970962052; 0.970962052; 0.970962052; 0.970962052;

0.970962052; 0.493270648; 0.493270648; 0.493270648;

0.493270648; 0.493270648; 0.493270648; 0.493270648;

0.493270648; 0.493270648; 0.493270648; 0.801578777;

0.801578777; 0.801578777; 0.801578777; 0.801578777;

0.801578777; 0.801578777; 0.801578777; 0.801578777;

0.794818776; 0.794818776; 0.794818776; 0.794818776;

0.794818776; 0.794818776; 0.794818776; 0.378519258;

0.378519258; 0.378519258; 0.378519258; 0.378519258;

0.378519258; 0.378519162];

end VFA_in;

record COD_s_nVFA_in

// soluble non-vfa COD

extends AMOCO.Dataset2020.TimeSeries;

parameter Real DataFrame_COD_s_nVFA_in[:, 2] =

[TimeSeries, Data_COD_s_nVFA_in];

parameter Real Data_COD_s_nVFA_in[:, 1] = [5.714927896;

5.714927896; 5.714927896; 5.714927896; 5.714927896;

6.310290225; 6.310290225; 6.310290225; 6.310290225;

6.905652553; 6.905652553; 6.905652553; 6.519592174;

6.519592174; 6.519592174; 6.133531795; 6.133531795;

6.133531795; 6.40556396; 6.40556396; 6.40556396; 6.40556396;

6.40556396; 6.40556396; 6.77966517; 6.77966517; 6.77966517;

6.77966517; 6.77966517; 6.854186607; 6.854186607; 6.854186607;

6.920909952; 7.254526677; 7.254526677; 7.254526677;

7.201248054; 7.201248054; 7.201248054; 7.011069198;

7.011069198; 7.011069198; 7.011069198; 6.956672382;

6.956672382; 6.977144287; 6.977144287; 7.077789655;

7.077789655; 7.077789655; 7.077789655; 7.077789655;

7.077789655; 7.64138326; 7.64138326; 7.64138326; 7.64138326;

7.64138326; 7.64138326; 7.64138326; 7.64138326; 7.64138326;

7.64138326; 7.64138326; 7.64138326; 7.64138326; 7.64138326;

7.64138326; 7.019543548; 7.019543548; 7.019543548;

7.019543548; 7.019543548; 7.019543548; 7.019543548;

7.019543548; 7.019543548; 7.019543548; 7.019543548;

6.397703835; 6.397703835; 6.397703835; 6.397703835;

6.397703835; 6.397703835; 6.397703835; 6.397703835;

6.680601214; 6.680601214; 6.680601214; 6.680601214;

6.680601214; 6.680601214; 6.680601214; 6.963498594;

6.963498594; 6.963498594; 6.963498594; 6.579856813;

6.579856813; 6.579856813; 6.579856813; 6.579856813;

8.52095035; 8.52095035; 8.52095035; 8.52095035; 8.52095035;

8.52095035; 4.328835188; 4.328835188; 4.328835188;

4.328835188; 4.328835188; 4.328835188; 4.328835188;

4.328835188; 4.328835188; 4.328835188; 7.03447982; 7.03447982;

7.03447982; 7.03447982; 7.03447982; 7.03447982; 7.03447982;

7.03447982; 7.03447982; 6.975155531; 6.975155531; 6.975155531;

6.975155531; 6.975155531; 6.975155531; 6.975155531;

3.321802124; 3.321802124; 3.321802124; 3.321802124;

3.321802124; 3.321802124; 3.321801282];

Appendix B

120

end COD_s_nVFA_in;

end Dataset2020;

model AMOCO_vanilla

// this model contains the invocation of parameters and

wariables which are needed

//for the equation execution.

//THis model does not work standalone because the input

concentrations for the COD,

//pH control and hydrogen injection are defined in

ControlRoom.

//> 150 d = 12960000 s

//extends QuiescentModelWithInheritance(gamma=0.3,

delta=0.01);

//inclusione parameteri e variabililo

extends AMOCO.Parameters.BioChemicalParameters;

extends AMOCO.Parameters.PhysioChemicalParameters;

extends AMOCO.Parameters.ReactorParameters;

extends AMOCO.Parameters.PartiotioningParameters;

extends AMOCO.Parameters.CODFractioningParameters;

//cambiare fra dataset rosito e dataset 2020 cambiare

anche le condizioni iniziali(se si cambia il dataset di

parametri

extends AMOCO.Parameters.InputParameters ;

extends AMOCO.Parameters.InhibitionParameters;

//"cambiare fra dataset rosito (InitialConditions) e

dataset 2020 - InitialConditions2020"

//extends AMOCO.Parameters.InitialConditions ;

//"invocare C.I.rosito/Dataset2020 nel record

StateVariables"

extends AMOCO.Variables.StateVariables ;

extends AMOCO.Variables.AlgebraicVariables;

extends AMOCO.Variables.RateVariables;

extends AMOCO.Variables.InhibitionVariables;

extends AMOCO.Variables.GasVariables;

// inclusione equazioni

/*parameteri reattore*/

Real h2_GORG_in(unit = "kmol/s", min = 0);

Real VALV_h2 ;

Real X_T "total solids ";

Real z_1;

Real z_2;

Real z_3;

equation

/*inhibition equations */

I_pH_1 = if noEvent(pH < pH_ul_1) then exp(-3 * ((pH -

pH_ul_1) / (pH_ul_1 - pH_ll_1)) ^ 2) else 1;

I_pH_2 = if noEvent(pH < pH_ul_2) then exp(-3 * ((pH -

pH_ul_2) / (pH_ul_2 - pH_ll_2)) ^ 2) else 1;

I_pH_3 = if noEvent(pH < pH_ul_3) then exp(-3 * ((pH -

pH_ul_3) / (pH_ul_3 - pH_ll_3)) ^ 2) else 1;

Appendix B

121

I_IC = 1 / (1 + K_I_IC / IC);

I_h2 = 1 / (1 + S_4 / K_I_h2);

/*rates*/

r_1 = k_h1 * X_01;

r_2 = mu_1_max / Y_x_1 * S_1 / (S_1 + K_s_1) * X_1 *

I_pH_1 * I_h2;

r_3 = mu_2_max / Y_x_2 * S_2 / (S_2 + K_s_2) * X_2 *

I_pH_2;

r_4 = mu_3_max / Y_x_3 * S_4 / (S_4 + K_s_4) * X_3 *

I_pH_3 * I_IC;

r_5 = k_d_1 * X_1;

r_6 = k_d_2 * X_2;

r_7 = k_d_3 * X_3;

r_8 = if noEvent(S_4 > 0) then k_la_h2 * (S_4 / 16 -

k_h_h2 * p_h2) else 0;

r_9 = if noEvent(S_5 > 0) then k_la * (S_5 / 64 - k_h_ch4

* p_ch4) else 0;

r_10 = if noEvent(S_3 > 0) then k_la * (S_3 - k_h_co2 *

p_co2) else 0;

r_11 = k_la_h2 * (k_h_h2 * p_gas - S_4 / 16) * (-VALV_h2);

r_12 = k_ab_co2 * (S_hco3 * (k_a_co2 + S_h) - k_a_co2 *

IC);

r_13 = k_ab_IN * (S_nh3 * (k_a_IN + S_h) - k_a_IN * IN);

/* differential equations */

der(X_01) = Q / V * (X_01_in_Signal - X_01) - r_1 + r_5 +

r_6 + r_7;

der(X_nd) = Q / V * (-X_nd) + k_nd * r_1;

der(X_1) = Q / V * (X_1_in - X_1) + Y_x_1 * r_2 - r_5;

der(X_2) = Q / V * (X_2_in - X_2) + Y_x_2 * r_3 - r_6;

der(X_3) = Q / V * (X_3_in - X_3) + Y_x_3 * r_4 - r_7;

der(S_1) = Q / V * (S_1_in_Signal - S_1) + (1 - k_nd) *

r_1 - r_2;

der(S_2) = Q / V * (S_2_in_Signal - S_2) + (1 - Y_x_1) *

f_vfa * r_2 - r_3;

der(S_4) = Q / V * (-S_4) + (1 - Y_x_1) * f_h2 * r_2 - r_4

- r_8 * 16 + (1 - f_bypass) * r_11 * 16;

der(S_5) = Q / V * (-S_5) + (1 - Y_x_2) * r_3 + (1 -

Y_x_3) * r_4 - r_9 * 64;

der(IC) = Q / V * (IC_in - IC) + r_1 * (C_x01 - C_s1) +

r_2 * (C_s1 - (Y_x_1 * C_bac + (1 - Y_x_1) * C_s2 * f_vfa)) +

r_3 * (C_s2 - (Y_x_2 * C_bac + (1 - Y_x_2) * C_ch4)) + r_4 *

(-(Y_x_3 * C_bac + (1 - Y_x_3) * C_ch4)) + (r_5 + r_6 + r_7) *

(C_bac - C_x01) - r_10;

der(S_hco3) = (-r_12) + Q / V * S_hco3_in_Signal;

der(S_nh3) = -r_13;

der(IN) = Q / V * (IN_in - IN) + r_1 * (N_x01 - N_s1) +

r_2 * (N_s1 - Y_x_1 * N_bac) + r_3 * (-Y_x_2 * N_bac) + r_4 *

(-Y_x_3 * N_bac) + (r_5 + r_6 + r_7) * (N_bac - N_x01);

der(S_cat) = Q / V * (S_cat_in - S_cat) - u_pH;

der(S_an) = Q / V * (S_an_in - S_an);

Appendix B

122

der(S_gas_h2) = (-q_gas * S_gas_h2 / V_gas) + V / V_gas *

r_8 + f_bypass * r_11 * 16;

/*-r_11*/

der(S_gas_ch4) = (-q_gas * S_gas_ch4 / V_gas) + V / V_gas

* r_9;

der(S_gas_co2) = (-q_gas * S_gas_co2 / V_gas) + V / V_gas

* r_10;

/* other equations*/

S_3 = IC - S_hco3;

S_nh4 = IN - S_nh3;

// idrogeno gorgogliato in ingresso (kgCOD/s)

h2_GORG_in = r_11 * V;

//alcalinità

Alk = (IC + k_w / S_h + S_2_ion + S_nh4 - S_h) * 1000 *

50;

/* ionic balance */

S_h = 10 ^ (-pH);

S_2_ion = k_a_ac * S_2 / 64 / (k_a_ac + S_h);

THETA = S_cat + S_nh4 - S_hco3 - S_an - S_2_ion / 64;

S_h = (-THETA / 2) + sqrt(THETA ^ 2 + 4 * k_w) / 2;

/* gas equations */

p_h2 = S_gas_h2 * 1000 * R * T;

p_ch4 = S_gas_ch4 * 1000 * R * T;

p_co2 = S_gas_co2 * 1000 * R * T;

p_tot = p_h2 + p_ch4 + p_co2;

q_gas = R * T / (p_gas - p_h2o) * V * 1000 * (r_8 + r_9 +

r_10);

x_h2 = p_h2 / p_tot;

x_ch4 = p_ch4 / p_tot;

x_co2 = p_co2 / p_tot;

//particolato toatale

X_T = X_01 + X_nd + X_1 + X_2 + X_3;

z_1 = (1 - Y_x_1) * f_vfa * r_2;

z_2 = - r_3;

der(z_3) = z_1 - z_2;

end AMOCO_vanilla;

model ControlRoom

// this is the main model which should be executed

// the

//20d = 1728000 s

//40d = 3456000 s

//140d = 12096000 s

//160d = 13824000 S

//180d = 15552000 s

extends AMOCO.AMOCO_vanilla;

extends AMOCO.DataSim;

//

//parameter definition

parameter Real S_hco3_Signal = 0 "for alkalinity control -

never implemented";

Appendix B

123

parameter Real pH_lim(unit = "Speece") = 8.2;

Real SludgeSignal "variable used to build an input signal

for COD";

Real VALV_h2 "variable regulated by a PI controller to

acheve the target level of hydrogen injection";

Real u_pH;

parameter Real co2_out_eql(unit = "kmol/s") = 6.23083E-10

"new:6.23083E-10 old:4.37635E-10";

//parameter Real X

Real h2_GORG_control(unit = "kmol/s");

Real pH_control;

Real eta;

Real X_bac(unit = "gCOD/l");

//

/* concentrations input definition - sludge signal can be

either user defined

or set equal to the dataset particualte COD input*/

AMOCO.Types.Conc_COD X_01_in_Signal =

/*COD_x_in_data;*/SludgeSignal*X_01_in;

AMOCO.Types.Conc_COD S_1_in_Signal = SludgeSignal*S_1_in

;

AMOCO.Types.Conc_COD S_2_in_Signal = SludgeSignal*S_2_in

;

AMOCO.Types.M S_hco3_in_Signal = S_hco3_start *

S_hco3_Signal * 1e3;

Real COD_T_in;

//

// hydrogen injection, ph, PI block call

AMOCO.Blocks.StepSignal GorgStep1;

AMOCO.Blocks.StepSignal GorgStep2;

AMOCO.Blocks.StepSignal GorgStep3;

AMOCO.Blocks.StepSignal GorgStep4;

AMOCO.Blocks.StepSignal GorgStep5;

AMOCO.Blocks.StepSignal GorgStep6;

AMOCO.Blocks.StepSignal GorgStep7;

AMOCO.Blocks.StepSignal pHStep1;

AMOCO.Blocks.StepSignal pHStep2;

AMOCO.Blocks.StepSignal pHStep3;

AMOCO.Blocks.ControllerPI PI;

AMOCO.Blocks.ControllerPI PI_pH;

//

//sludge modulation step block invocation

AMOCO.Blocks.StepSignal SludgeStep1;

AMOCO.Blocks.StepSignal SludgeStep2;

AMOCO.Blocks.StepSignal SludgeStep3;

AMOCO.Blocks.StepSignal SludgeStep4;

//

// (USER INPUT REQUIRED)reactor control:

// by setting true or false the ph, hydrogen injection and

sludge

// input control is enabled

Appendix B

124

parameter Boolean CONTROL_pH = false;

parameter Boolean CONTROL_GORG = false;

parameter Boolean CONTROL_SLUDGE = false ;

initial equation

PI_pH.y = 0;

equation

//

//ph control signal construction

pHStep1.Offset = pH;

pHStep1.StepTime = 86400 * 20;

pHStep1.Step = 7.3;

connect(pHStep2.Offset, pHStep1.y);

pHStep2.StepTime = 86400 * 40;

pHStep2.Step = 7.3;

connect(pHStep3.Offset, pHStep2.y);

pHStep3.StepTime = 86400 * 100;

pHStep3.Step = 7.3;

connect(pH_control, pHStep3.y);

//

//ph pi controller implementation

PI_pH.K_p = 1e-5;

PI_pH.K_i = 1e-10;

connect(PI_pH.u, pH);

connect(PI_pH.c, pH_control) "pH_control";

if CONTROL_pH == true then

connect(PI_pH.y, u_pH);

else

u_pH = 0;

end if;

//

//Pi control connection for hydrogen injection or GORG valv

PI.K_p = 1000e1 "proportional gain";

PI.K_i = 200e1 "integral time";

connect(PI.u, h2_GORG_in);

connect(PI.c, h2_GORG_control);

if CONTROL_GORG == true then

connect(PI.y, VALV_h2);

else

VALV_h2 = 0;

end if;

//

// hydrogen injection target signal definition

GorgStep1.Offset = 0;

GorgStep1.StepTime = 86400 * 21;

GorgStep1.Step = co2_out_eql * 0.76 "0.76";

connect(GorgStep1.y, GorgStep2.Offset);

GorgStep2.StepTime = 86400 * 42;

GorgStep2.Step = co2_out_eql * 1.01"1.01";

connect(GorgStep2.y, GorgStep3.Offset);

GorgStep3.StepTime = 86400 * 50;

GorgStep3.Step = co2_out_eql * 1.9"1.9";

Appendix B

125

connect(GorgStep3.y, GorgStep4.Offset);

GorgStep4.StepTime = 86400 * 59;

GorgStep4.Step = co2_out_eql * 2.6"2.6";

connect(GorgStep4.y, GorgStep5.Offset);

GorgStep5.StepTime = 86400 * 68;

GorgStep5.Step = co2_out_eql * 4.36 "4.36";

connect(GorgStep5.y, GorgStep6.Offset);

GorgStep6.StepTime = 86400 * 95;

GorgStep6.Step = co2_out_eql * 6.01 "6.01";

connect(GorgStep6.y, GorgStep7.Offset);

GorgStep7.StepTime = 86400 * 108;

GorgStep7.Step = co2_out_eql * 7 "7";

connect(GorgStep7.y, h2_GORG_control);

//

// Sludge input Modulation

// set each step to 1 to disable all modifications to the

input signal

SludgeStep1.Offset = 1;

SludgeStep1.StepTime = 86400 * 20;

SludgeStep1.Step = (1);

connect(SludgeStep1.y, SludgeStep2.Offset);

SludgeStep2.StepTime = 86400 * 40;

SludgeStep2.Step = (1);

connect(SludgeStep2.y, SludgeStep3.Offset);

SludgeStep3.StepTime = 86400 * 60;

SludgeStep3.Step = (1);

connect(SludgeStep3.y, SludgeStep4.Offset);

SludgeStep4.StepTime = 86400 * 80;

SludgeStep4.Step = (1);

connect(SludgeStep4.y, SludgeSignal);

// total particulate cod input

COD_T_in = X_01_in_Signal + S_1_in_Signal + S_2_in_Signal;

// dummy code string - to be removed

eta = 1 - (X_01_in_Signal + S_1_in_Signal + S_2_in_Signal

+ X_1_in + X_2_in + X_3_in)/(X_01 + X_nd + X_1 + X_2 + X_3 +

S_1 + S_2 + S_4 + S_5 );

X_bac = X_1 + X_2 + X_3;

end ControlRoom;

package Blocks

//this package contains the blocks implemented in the code

block StepSignal

// block useful to create step signals

Real Offset "starting value";

//Real TimeInput;

Real Step "step value";

Real StepTime "time at which the step occurs";

Real y;

equation

y = if time < StepTime then Offset else Step;

Appendix B

126

end StepSignal;

block ControllerPI

//block used to create PI-controlled signals

Real u "input signal";

Real c "control signal";

Real y "output signal";

Real K_p "proportional gain";

Real K_i "integral time";

protected

Real e "error";

equation

e = u - c;

der(y) = K_p * der(e) + K_i * e;

end ControllerPI;

end Blocks;

model DataSim

//modello che prende le dataseries e le plotta

// invocando tutti i record presenti nel package del

dataset

//è utile per confrontare i dati e definire entitiò

dell'errore

// i valori di input possono essere passati al modello se

non si

// vuole usare il valore medio

//40d = 3456000 s

//sim time = 12182400 (160 d)

// sim time = 12096000 (140d)

extends AMOCO.Dataset2020.Alk;

extends AMOCO.Dataset2020.VFA;

extends AMOCO.Dataset2020.CodSol;

extends AMOCO.Dataset2020.Q_gas;

extends AMOCO.Dataset2020.x_co2;

extends AMOCO.Dataset2020.x_ch4;

extends AMOCO.Dataset2020.x_h2;

extends AMOCO.Dataset2020.COD_x_out;

extends AMOCO.Dataset2020.pH;

//invocazione input

extends AMOCO.Dataset2020.COD_x_in;

extends AMOCO.Dataset2020.VFA_in;

extends AMOCO.Dataset2020.COD_s_nVFA_in;

//definizione varibili delle tavole

// timetable crea una variabile tempo dipendente

// a partire da una matrice

//variabili di output

Modelica.Blocks.Sources.TimeTable AlkSim(table =

DataFrameAlk);

Modelica.Blocks.Sources.TimeTable VFASim(table =

DataFrameVFA);

Appendix B

127

Modelica.Blocks.Sources.TimeTable CODSim(table =

DataFrameCodSol);

Modelica.Blocks.Sources.TimeTable Q_gasSim(table =

DataFrameQ_gas);

Modelica.Blocks.Sources.TimeTable x_co2_Sim(table =

DataFrame_x_co2);

Modelica.Blocks.Sources.TimeTable x_ch4_Sim(table =

DataFrame_x_ch4);

Modelica.Blocks.Sources.TimeTable x_h2_Sim(table =

DataFrame_x_h2);

Modelica.Blocks.Sources.TimeTable COD_x_out_Sim(table =

DataFrame_COD_x_out);

Modelica.Blocks.Sources.TimeTable pHSim(table =

DataFramepH);

//input

Modelica.Blocks.Sources.TimeTable COD_x_in_Sim(table =

DataFrame_COD_x_in);

Modelica.Blocks.Sources.TimeTable VFA_in_Sim(table =

DataFrame_VFA_in);

Modelica.Blocks.Sources.TimeTable COD_s_nVFA_in_Sim(table

= DataFrame_COD_s_nVFA_in);

//dichiarazione di variabili per le serie temporali -

dichiarate per

// semplificare la ricerca nel plot window

Real Alk_data(unit = "mgCaCO3/m3");

Real VFA_data(unit = "mgCOD/l");

Real COD_s_tot_data(unit = "mgCOD/l");

Real COD_x_out_data(unit = "gCOD/l");

Real Q_gas_data(unit = "ml/d");

Real x_co2_data;

Real x_ch4_data;

Real x_h2_data;

Real pH_data;

Real COD_x_in_data;

Real VFA_in_data(unit = "kgCOD/m3");

Real COD_s_nVFA_in_data(unit = "kgCOD/m3");

//

// derived variables

Real q_h2_data(unit ="ml/d");

Real q_ch4_data(unit ="ml/d");

Real q_co2_data(unit ="ml/d");

Real S_1_data(unit="mgCOD/l");

Real S_2_data(unit="mgCOD/l");

equation

// connessione fra le tabelle e le variabili

connect(Alk_data, AlkSim.y);

connect(VFA_data, VFASim.y);

connect(COD_s_tot_data, CODSim.y);

connect(Q_gas_data, Q_gasSim.y);

connect(x_co2_data, x_co2_Sim.y);

Appendix B

128

connect(x_ch4_data, x_ch4_Sim.y);

connect(x_h2_data, x_h2_Sim.y);

connect(COD_x_out_data, COD_x_out_Sim.y);

connect(pH_data, pHSim.y);

connect(COD_x_in_data, COD_x_in_Sim.y);

connect(VFA_in_data, VFA_in_Sim.y);

connect(COD_s_nVFA_in_data, COD_s_nVFA_in_Sim.y);

//

//gas flow - creazione di serie temporali per il

frazionamento

//gassoso e del COD solubile

q_h2_data = Q_gas_data*x_h2_data/100;

q_ch4_data = Q_gas_data*x_ch4_data/100;

q_co2_data = Q_gas_data*x_co2_data/100;

S_1_data = COD_s_tot_data - VFA_data;

S_2_data = VFA_data;

end DataSim;

model DataConversion

//questo modello:

// (1) prende in considerazione le variablil esistenti e le

converte da

// un formato strettamente SI a uno piu leggibile

// (2) calcola il bilancio di masssa

// 40d = 3456000 s

// sim time = 13824000 (160 d)

// sim time = 12096000 (140d)

extends AMOCO.ControlRoom;

Real Alk_Model(unit = "mgCaCO3/m3");

Real VFA_Model(unit = "mgHAC/l");

Real COD_sol_bio_out_Model(unit = "mgCOD/l");

Real Q_gas_out_Model(unit = "ml/d");

Real q_h2_model(unit = "ml/d");

Real q_ch4_model(unit = "ml/d");

Real q_co2_model(unit = "ml/d");

Real q_h2_in_Signal(unit="ml/d");

Real q_h2_in(unit="ml/d");

Real x_co2_out_Model;

Real x_ch4_out_Model;

Real pH_out_Model;

// mgCOD/l conversion & mmol/l &

AMOCO.Types.Conc_COD_model X_01_model;

AMOCO.Types.Conc_COD_model X_nd_model;

AMOCO.Types.Conc_COD_model X_1_model;

AMOCO.Types.Conc_COD_model X_2_model;

AMOCO.Types.Conc_COD_model X_3_model;

AMOCO.Types.Conc_COD_model S_1_model;

AMOCO.Types.Conc_COD_model S_2_model;

AMOCO.Types.Conc_COD_model S_4_model;

AMOCO.Types.Conc_COD_model S_5_model;

AMOCO.Types.M_model IC_model;

Appendix B

129

AMOCO.Types.M_model IN_model;

AMOCO.Types.M_model S_hco3_model;

AMOCO.Types.M_model S_co2_model;

AMOCO.Types.M_model S_nh3_model;

AMOCO.Types.M_model S_nh4_model;

AMOCO.Types.M_model S_cat_model;

AMOCO.Types.M_model S_an_model;

AMOCO.Types.M_model S_gas_h2_model;

AMOCO.Types.M_model S_gas_ch4_model;

AMOCO.Types.M_model S_gas_co2_model;

//

//COD mass balance variable declaration

AMOCO.Types.MassFlowDay m_dot_X_01_in;

AMOCO.Types.MassFlowDay m_dot_X_1_in;

AMOCO.Types.MassFlowDay m_dot_X_2_in;

AMOCO.Types.MassFlowDay m_dot_X_3_in;

AMOCO.Types.MassFlowDay m_dot_S_1_in;

AMOCO.Types.MassFlowDay m_dot_S_2_in;

AMOCO.Types.MassFlowDay m_dot_S_gas_h2_GORG;

AMOCO.Types.MassFlowDay m_dot_X_01;

AMOCO.Types.MassFlowDay m_dot_X_nd;

AMOCO.Types.MassFlowDay m_dot_X_1;

AMOCO.Types.MassFlowDay m_dot_X_2;

AMOCO.Types.MassFlowDay m_dot_X_3;

AMOCO.Types.MassFlowDay m_dot_S_1;

AMOCO.Types.MassFlowDay m_dot_S_2;

AMOCO.Types.MassFlowDay m_dot_S_4;

AMOCO.Types.MassFlowDay m_dot_S_5;

AMOCO.Types.MassFlowDay m_dot_S_gas_h2;

AMOCO.Types.MassFlowDay m_dot_S_gas_ch4;

AMOCO.Types.MassFlowDay COD_in;

AMOCO.Types.MassFlowDay COD_out;

Real COD_balance(unit ="kgCOD/d");

//frazionamento gassoso come rapporto fra substrati

Real x_ch4_2;

Real x_co2_2;

Real x_h2_2;

//

//olr

Real OLR(unit="kgCOD/d/m3");

equation

Alk_Model = Alk;

VFA_Model = S_2 * 1000;

COD_sol_bio_out_Model = S_1 * 1000;

Q_gas_out_Model = q_gas * 86400 * 1e6;

q_h2_model = Q_gas_out_Model*x_h2;

q_ch4_model = Q_gas_out_Model*x_ch4;

q_co2_model = Q_gas_out_Model*x_co2;

q_h2_in_Signal = h2_GORG_control * R*T/p_gas *86400*1e9;

q_h2_in = h2_GORG_in * R*T/p_gas *86400*1e9;

x_co2_out_Model = x_co2 * 100;

Appendix B

130

x_ch4_out_Model = x_ch4 * 100;

pH_out_Model = pH;

//

//concentration in terms of mmol/l mg/l

X_01_model = X_01 * 1000 ;

X_nd_model = X_nd * 1000 ;

X_1_model = X_1 * 1000 ;

X_2_model = X_2 * 1000 ;

X_3_model = X_3 * 1000 ;

S_1_model = S_1 * 1000 ;

S_2_model = S_2 * 1000 ;

S_4_model = S_4 * 1000 ;

S_5_model = S_5 * 1000 ;

IC_model = IC * 1000 ;

IN_model = IN * 1000 ;

S_hco3_model= S_hco3 * 1000 ;

S_co2_model = S_3 * 1000 ;

S_nh3_model = S_nh3 * 1000 ;

S_nh4_model = S_nh4 * 1000 ;

S_cat_model = S_cat * 1000 ;

S_an_model = S_an * 1000 ;

S_gas_h2_model = S_gas_h2 * 1000 ;

S_gas_ch4_model = S_gas_ch4 * 1000 ;

S_gas_co2_model = S_gas_co2 * 1000 ;

//

// daily mass flows in grams

////input

m_dot_X_01_in = X_01_in_Signal * Q * 1000 * 86400;

m_dot_X_1_in = X_1_in * Q * 1000 * 86400;

m_dot_X_2_in = X_2_in * Q * 1000 * 86400;

m_dot_X_3_in = X_3_in * Q * 1000 * 86400;

m_dot_S_1_in = S_1_in_Signal * Q * 1000 * 86400;

m_dot_S_2_in = S_2_in_Signal * Q * 1000 * 86400;

m_dot_S_gas_h2_GORG = h2_GORG_in * 16 * 1000 * 86400;

//

////output

m_dot_X_01 = X_01 * Q * 1000 * 86400;

m_dot_X_nd = X_nd * Q * 1000 * 86400;

m_dot_X_1 = X_1 * Q * 1000 * 86400;

m_dot_X_2 = X_2 * Q * 1000 * 86400;

m_dot_X_3 = X_3 * Q * 1000 * 86400;

m_dot_S_1 = S_1 * Q * 1000 * 86400;

m_dot_S_2 = S_2 * Q * 1000 * 86400;

m_dot_S_4 = S_4 * Q * 1000 * 86400;

m_dot_S_5 = S_5 * Q * 1000 * 86400;

m_dot_S_gas_h2 = S_gas_h2 * q_gas * 16 * 1000 * 86400;

m_dot_S_gas_ch4 = S_gas_ch4 * q_gas * 64 * 1000 * 86400;

//

//mass balance

COD_in = m_dot_X_01_in

+ m_dot_X_1_in

Appendix B

131

+ m_dot_X_2_in

+ m_dot_X_3_in

+ m_dot_S_1_in

+ m_dot_S_2_in

+ m_dot_S_gas_h2_GORG;

COD_out = m_dot_X_01

+ m_dot_X_nd

+ m_dot_X_1

+ m_dot_X_2

+ m_dot_X_3

+ m_dot_S_1

+ m_dot_S_2

+ m_dot_S_4

+ m_dot_S_5

+ m_dot_S_gas_h2

+ m_dot_S_gas_ch4;

COD_balance = COD_in - COD_out;

x_ch4_2 = S_gas_ch4 / (S_gas_ch4 + S_gas_co2 + S_gas_h2);

x_co2_2 = S_gas_co2 / (S_gas_ch4 + S_gas_co2 + S_gas_h2);

x_h2_2 = S_gas_h2 / (S_gas_ch4 + S_gas_co2 + S_gas_h2);

//

//OLR

OLR = COD_in/V;

end DataConversion;

model Carbon_Nitrogen_Balance

// questo modello calcola i bilanci di massa di carbonio e

azoto

//

//object import

extends AMOCO.Parameters.PartiotioningParameters;

extends AMOCO.ControlRoom;

//

//parameters

parameter Real C_xnd = 0.05 "default = 0.02786" ;

parameter Real N_I = 0.06/14 ;

//

//carbon input variables

Real C_dot_X_01_in(unit="molC/d") ;

Real C_dot_X_bac_in(unit="molC/d") ;

Real C_dot_S_1_in(unit="molC/d") ;

Real C_dot_S_2_in(unit="molC/d") ;

Real C_dot_IC_in(unit="molC/d") ;

//

//Nitrogen Input Variables

Real N_dot_X_01_in(unit="molN/d") ;

Real N_dot_X_bac_in(unit="molN/d") ;

Real N_dot_S_1_in(unit="molN/d") ;

Real N_dot_IN_in(unit="molN/d") ;

//

//carbon output variables

Appendix B

132

Real C_dot_X_01(unit="molC/d") ;

Real C_dot_X_nd(unit="molC/d") ;

Real C_dot_X_bac(unit="molC/d") ;

Real C_dot_S_1(unit="molC/d") ;

Real C_dot_S_2(unit="molC/d") ;

Real C_dot_S_5(unit="molC/d") ;

Real C_dot_IC(unit="molC/d") ;

Real C_dot_S_gas_ch4(unit="molC/d") ;

Real C_dot_S_gas_co2(unit="molC/d") ;

//

//Nitrogen Output Variables

Real N_dot_X_01(unit="molN/d") ;

Real N_dot_X_nd(unit="molN/d") ;

Real N_dot_X_bac(unit="molN/d") ;

Real N_dot_S_1(unit="molN/d") ;

Real N_dot_IN(unit="molN/d") ;

//

// Carbon Balance variables

Real C_dot_in(unit="molC/d");

Real C_dot_out(unit="molC/d");

Real C_dot_net(unit="molC/d");

//

// Nitrogen Balance Equations

Real N_dot_in(unit="molN/d");

Real N_dot_out(unit="molN/d");

Real N_dot_net(unit="molN/d");

//

equation

//

// carbon input equations

C_dot_X_01_in = (Q*86400*1000)*X_01_in*C_x01 ;

C_dot_X_bac_in = (Q*86400*1000)*(X_2_in + X_2_in +

X_3_in)*C_bac ;

C_dot_S_1_in = (Q*86400*1000)*S_1_in*C_s1 ;

C_dot_S_2_in = (Q*86400*1000)*S_2_in*C_s2 ;

C_dot_IC_in = (Q*86400*1000)*IC_in ;

//

// carbon output equations

C_dot_X_01 = (Q*86400*1000)*X_01*C_x01 ;

C_dot_X_nd = (Q*86400*1000)*X_nd*C_xnd ;

C_dot_X_bac = (Q*86400*1000)*(X_1 + X_2 + X_3)*C_bac ;

C_dot_S_1 = (Q*86400*1000)*S_1*C_s1 ;

C_dot_S_2 = (Q*86400*1000)*S_2*C_s2 ;

C_dot_S_5 = (Q*86400*1000)*S_5*C_ch4 ;

C_dot_IC = (Q*86400*1000)*IC ;

C_dot_S_gas_ch4 = (q_gas*86400*1000)*S_gas_ch4 ;

C_dot_S_gas_co2 = (q_gas*86400*1000)*S_gas_co2 ;

//

//balance equations

C_dot_in = C_dot_X_01_in + C_dot_X_bac_in + C_dot_S_1_in +

C_dot_S_2_in + C_dot_IC_in ;

Appendix B

133

C_dot_out = C_dot_X_01 + C_dot_X_nd + C_dot_X_bac +

C_dot_S_1 + C_dot_S_2 + C_dot_S_5 + C_dot_IC + C_dot_S_gas_ch4

+ C_dot_S_gas_co2 ;

C_dot_net = C_dot_in - C_dot_out ;

//

// Nitrogen input equations

N_dot_X_01_in = (Q*86400*1000)*X_01_in*N_x01 ;

N_dot_X_bac_in = (Q*86400*1000)*(X_1_in + X_2_in +

X_3_in)*N_bac ;

N_dot_S_1_in = (Q*86400*1000)*S_1_in*N_s1 ;

N_dot_IN_in = (Q*86400*1000)*IN_in ;

//

// Nitrogen Output equations

N_dot_X_01 = (Q*86400*1000)*X_01*N_x01 ;

N_dot_X_nd = (Q*86400*1000)*X_nd*N_I ;

N_dot_X_bac = (Q*86400*1000)*(X_1 + X_2 + X_3)*N_bac ;

N_dot_S_1 = (Q*86400*1000)*S_1 ;

N_dot_IN = (Q*86400*1000)*IN ;

//

// Nitrogen Balance Equations

N_dot_in = N_dot_X_01_in + N_dot_X_bac_in + N_dot_S_1_in +

N_dot_IN_in ;

N_dot_out = N_dot_X_01 + N_dot_X_nd + N_dot_X_bac +

N_dot_S_1 + N_dot_IN ;

N_dot_net = N_dot_in - N_dot_out ;

//

end Carbon_Nitrogen_Balance;

end AMOCO;

APPENDIX C

Matlab Code

C.1 LFTfun Function

% lftfun function definition

% (activ) SECTION - user input required

% (PASSIVE) SECTION - user input not required

%

% (ACTIVE) define parameters

function [lftfun] = AMOCO_2_18_lft(R,...

T,...

V,...

V_gas,...

Q,...

p_gas,...

p_h2o,...

mu_1_max,...

mu_2_max,...

mu_3_max,...

k_d_1,...

k_d_2,...

k_d_3,...

K_s_1,...

K_s_2,...

K_s_4,...

Y_x_1,...

Y_x_2,...

Y_x_3,...

k_h1,...

k_nd,...

C_bac,...

C_x01,...

C_s1,...

C_s2,...

C_ch4,...

Appendix C

135

N_bac,...

N_x01,...

N_s1,...

k_ab_co2,...

k_a_co2,...

k_ab_IN,...

k_a_IN,...

k_w,...

k_h_h2,...

k_h_ch4,...

k_h_co2,...

f_vfa_lim,...

k_la_h2_lim,...

k_la_lim)

%

%% (ACTIVE) define u, x, y, om numerosity

u_count = 10;

x_count = 18;

y_count = 9;

om_count = 18;

%

%% (ACTIVE) delta matrix definition

DeltaSym = {

'parName' 'indStartDiag' 'indStopDiag' 'LowerBound'

'HigherBound' 'toIdentify', 'lb', 'ub';

'f_vfa' 1 1 f_vfa_lim(1) f_vfa_lim(2) 1 -1 +1;

'k_la_h2' 2 5 k_la_h2_lim(1) k_la_h2_lim(2) 1 -1 +1;

'k_la' 6 13 k_la_lim(1) k_la_lim(2) 1 -1 +1;

};

delta_count = max(cell2mat(DeltaSym(2:end,3)));

%

%% (ACTIVE) THETA DEFINITION

om = sym('om', om_count);

theta_sym = [

(om(3)*om(6))/(K_s_1 + om(6))

(om(4)*om(7))/(K_s_2 + om(7))

(om(5)*om(8))/(K_s_4 + om(8))

om(12)*(om(7)/128 - om(11)/2 + om(12)/2 + om(13)/2 -

om(14)/2 + om(15)/2 + (4*k_w + (om(7)/64 - om(11) + om(12) +

om(13) - om(14) + om(15))^2)^(1/2)/2)

om(13)*(om(7)/128 - om(11)/2 + om(12)/2 + om(13)/2 -

om(14)/2 + om(15)/2 + (4*k_w + (om(7)/64 - om(11) + om(12) +

om(13) - om(14) + om(15))^2)^(1/2)/2)

om(16)*(om(8)/16 - R*T*k_h_h2*om(16))

om(17)*(om(8)/16 - R*T*k_h_h2*om(16))

om(18)*(om(8)/16 - R*T*k_h_h2*om(16))

om(16)*(om(9)/64 - R*T*k_h_ch4*om(17))

om(17)*(om(9)/64 - R*T*k_h_ch4*om(17))

om(18)*(om(9)/64 - R*T*k_h_ch4*om(17))

-om(16)*(om(12) - om(10) + R*T*k_h_co2*om(18))

-om(17)*(om(12) - om(10) + R*T*k_h_co2*om(18))

Appendix C

136

-om(18)*(om(12) - om(10) + R*T*k_h_co2*om(18))

(om(7)/128 - om(11)/2 + om(12)/2 + om(13)/2 - om(14)/2 +

om(15)/2 + (4*k_w + (om(7)/64 - om(11) + om(12) + om(13) -

om(14) + om(15))^2)^(1/2)/2)

1/((om(7)/128 - om(11)/2 + om(12)/2 + om(13)/2 - om(14)/2

+ om(15)/2 + (4*k_w + (om(7)/64 - om(11) + om(12) + om(13) -

om(14) + om(15))^2)^(1/2)/2))

];

theta_count = size(theta_sym,1);

%

%% DEFINITION OF THE DIMENSIONS OF LFT MATRICES (PASSIVE)

LTI = struct(...

'A', zeros(x_count,x_count),...

'B1', zeros(x_count,delta_count),...

'B2', zeros(x_count,theta_count),...

'B3', zeros(x_count,u_count),...

'C1', zeros(delta_count,x_count),...

'D11', zeros(delta_count,delta_count),...

'D12', zeros(delta_count,theta_count),...

'D13', zeros(delta_count,u_count),...

'C2', zeros(om_count,x_count),...

'D21', zeros(om_count,delta_count),...

'D22', zeros(om_count,theta_count),...

'D23', zeros(om_count,u_count),...

'C3', zeros(y_count,x_count),...

'D31', zeros(y_count,delta_count),...

'D32', zeros(y_count,theta_count),...

'D33', zeros(y_count,u_count));

%

%% (passive) MATRICES DEFINITION

% X_dot = A*x + B1*w + B2*XI + B3*u

syms x_ [1,x_count]

syms w_ [1,delta_count]

syms XI_ [1,theta_count]

syms u_ [1,u_count]

%% (ACTIVE) write LFT equations

%% nota bene: ciascuna equazione deve essere compresa in una

% sola riga

eqns = [

k_d_1*x_3 + k_d_2*x_4 + k_d_3*x_5 - k_h1*x_1 + (Q*u_1)/V -

(Q*x_1)/V

k_h1*k_nd*x_1 - (Q*x_2)/V

XI_1*mu_1_max - k_d_1*x_3 + (Q*u_2)/V - (Q*x_3)/V

XI_2*mu_2_max - k_d_2*x_4 + (Q*u_3)/V - (Q*x_4)/V

XI_3*mu_3_max - k_d_3*x_5 + (Q*u_4)/V - (Q*x_5)/V

k_h1*x_1 - k_h1*k_nd*x_1 + (Q*u_5)/V - (XI_1*mu_1_max)/Y_x_1 -

(Q*x_6)/V

Appendix C

137

(Q*u_6)/V - mu_1_max*w_1 - (XI_2*mu_2_max)/Y_x_2 - (Q*x_7)/V +

(mu_1_max*w_1)/Y_x_1

mu_1_max*w_1 - XI_1*mu_1_max - 16*w_2 + (XI_1*mu_1_max)/Y_x_1

- (XI_3*mu_3_max)/Y_x_3 - (Q*x_8)/V - (mu_1_max*w_1)/Y_x_1

(XI_2*mu_2_max)/Y_x_2 - XI_2*mu_2_max - XI_3*mu_3_max - 64*w_6

+ (XI_3*mu_3_max)/Y_x_3 - (Q*x_9)/V

C_ch4*XI_2*mu_2_max - C_bac*XI_1*mu_1_max -

C_bac*XI_2*mu_2_max - C_bac*XI_3*mu_3_max - w_7 +

C_ch4*XI_3*mu_3_max + C_bac*k_d_1*x_3 + C_bac*k_d_2*x_4 +

C_bac*k_d_3*x_5 - C_s1*k_h1*x_1 - C_x01*k_d_1*x_3 -

C_x01*k_d_2*x_4 - C_x01*k_d_3*x_5 + C_x01*k_h1*x_1 +

C_s2*mu_1_max*w_1 + (Q*u_7)/V - (Q*x_10)/V -

(C_ch4*XI_2*mu_2_max)/Y_x_2 - (C_ch4*XI_3*mu_3_max)/Y_x_3 +

(C_s1*XI_1*mu_1_max)/Y_x_1 + (C_s2*XI_2*mu_2_max)/Y_x_2 -

(C_s2*mu_1_max*w_1)/Y_x_1

N_bac*k_d_1*x_3 - N_bac*XI_2*mu_2_max - N_bac*XI_3*mu_3_max -

N_bac*XI_1*mu_1_max + N_bac*k_d_2*x_4 + N_bac*k_d_3*x_5 -

N_s1*k_h1*x_1 - N_x01*k_d_1*x_3 - N_x01*k_d_2*x_4 -

N_x01*k_d_3*x_5 + N_x01*k_h1*x_1 + (Q*u_8)/V - (Q*x_11)/V +

(N_s1*XI_1*mu_1_max)/Y_x_1

k_a_co2*k_ab_co2*x_10 - XI_4*k_ab_co2 - k_a_co2*k_ab_co2*x_12

k_a_IN*k_ab_IN*x_11 - XI_5*k_ab_IN - k_a_IN*k_ab_IN*x_13

(Q*u_9)/V - (Q*x_14)/V

(Q*u_10)/V - (Q*x_15)/V

(V*w_2)/V_gas - (1000*R*T*V*w_3)/(V_gas*p_gas - V_gas*p_h2o) -

(1000*R*T*V*w_8)/(V_gas*p_gas - V_gas*p_h2o) -

(1000*R*T*V*w_11)/(V_gas*p_gas - V_gas*p_h2o)

(V*w_6)/V_gas - (1000*R*T*V*w_4)/(V_gas*p_gas - V_gas*p_h2o) -

(1000*R*T*V*w_9)/(V_gas*p_gas - V_gas*p_h2o) -

(1000*R*T*V*w_12)/(V_gas*p_gas - V_gas*p_h2o)

(V*w_7)/V_gas - (1000*R*T*V*w_5)/(V_gas*p_gas - V_gas*p_h2o) -

(1000*R*T*V*w_10)/(V_gas*p_gas - V_gas*p_h2o) -

(1000*R*T*V*w_13)/(V_gas*p_gas - V_gas*p_h2o)

];

%

%% (PASSIVE) creation of the x_dot array

vars = [x w XI u];

[sym_E2M] = equationsToMatrix(eqns,vars);

E2M = double(sym_E2M);

%

% E2M matriceses dimensions definition

x_start=1;

x_stop=length(x);

w_start=length(x)+1;

w_stop=length(x)+length(w);

XI_start=length(x)+length(w)+1;

XI_stop= length(x)+length(w)+length(XI);

u_start=XI_stop+1;

%

% fractioning

A_E2M = E2M(:,x_start:x_stop);

Appendix C

138

B1_E2M = E2M(:,w_start:w_stop);

B2_E2M = E2M(:,XI_start:XI_stop);

B3_E2M = E2M(:,u_start:end);

%

%from E2M TI LTI

LTI.A=A_E2M;

LTI.B1=B1_E2M;

LTI.B2=B2_E2M;

LTI.B3=B3_E2M;

%

%% Z ARRAY CREATION

% Z = C1*x + D11*w + D12*XI + D13*u

%

%% (ACTIVE) z equations input

eqns_z = [

XI_1

x_8/16 - R*T*k_h_h2*x_16

XI_6

XI_7

XI_8

x_9/64 - R*T*k_h_ch4*x_17

x_10 - x_12 - R*T*k_h_co2*x_18

XI_9

XI_10

XI_11

XI_12

XI_13

XI_14

];

%

%% (PASSIVE) CREATION OF THE Z ARRAY FOR THE lft STRUCT

[sym_E2M_z] = equationsToMatrix(eqns_z,vars);

E2M_z = double(sym_E2M_z);

%

%fractioning(PASSIVE)

C1_E2M = E2M_z(:,x_start:x_stop);

D11_E2M = E2M_z(:,w_start:w_stop);

D12_E2M = E2M_z(:,XI_start:XI_stop);

D13_E2M = E2M_z(:,u_start:end);

%

%from E2M TI LTI (PASSIVE)

LTI.C1=C1_E2M;

LTI.D11=D11_E2M;

LTI.D12=D12_E2M;

LTI.D13=D13_E2M;

%

%% OM ARRAY CREATION

%OM = C2*x + D21*w + D22*XI + D23*u

%

%% (ACTIVE) OM equations input

eqns_om = [

Appendix C

139

x_1

x_2

x_3

x_4

x_5

x_6

x_7

x_8

x_9

x_10

x_11

x_12

x_13

x_14

x_15

x_16

x_17

x_18

];

%

%% (PASSIVE) CREATION OF THE OM ARRAY FOR THE lft STRUCT

[sym_E2M_om] = equationsToMatrix(eqns_om,vars);

E2M_om = double(sym_E2M_om);

%

%fractioning

C2_E2M = E2M_om(:,x_start:x_stop);

D21_E2M = E2M_om(:,w_start:w_stop);

D22_E2M = E2M_om(:,XI_start:XI_stop);

D23_E2M = E2M_om(:,u_start:end);

%

%from E2M TI LTI (PASSIVE)

LTI.C2=C2_E2M;

LTI.D21=D21_E2M;

LTI.D22=D22_E2M;

LTI.D23=D23_E2M;

%

%% Y ARRAY CREATION

% Y = C3*x + D31*w + D32*XI + D33*u

%

%% (ACTIVE) Y equations input

%

eqns_y = [

x_1 + x_2 + x_3 + x_4 + x_5 %solidi totale

x_6 %cod_s

x_7 %vfa

(x_7/64 - XI_15 + x_10 + x_11 - x_13 + k_w*XI_16)*1000*50

%alk

R * T / (p_gas - p_h2o) * V * 1000 * (w_2 + w_6 + w_7)

x_16

x_17

Appendix C

140

x_18

XI_15

];

%

%% (PASSIVE) CREATION OF THE Z ARRAY FOR THE lft STRUCT

[sym_E2M_y] = equationsToMatrix(eqns_y,vars);

E2M_y = double(sym_E2M_y);

%

%fractioning

C3_E2M = E2M_y(:,x_start:x_stop);

D31_E2M = E2M_y(:,w_start:w_stop);

D32_E2M = E2M_y(:,XI_start:XI_stop);

D33_E2M = E2M_y(:,u_start:end);

%

%from E2M TI LTI (PASSIVE)

LTI.C3=C3_E2M;

LTI.D31=D31_E2M;

LTI.D32=D32_E2M;

LTI.D33=D33_E2M;

%

%% LFTFUN DEFINITION (PASSIVE)

lftfun = struct(...

'LTI', LTI,...

'DeltaSym', {DeltaSym},...

'DeltaVal', zeros(delta_count,delta_count)...

);

%

% save additional data for reference purpose

lftfun.theta_sym = theta_sym;

lftfun.u_count = u_count;

lftfun.x_count = x_count;

lftfun.y_count = y_count;

lftfun.om_count = om_count;

lftfun.theta_count = theta_count;

lftfun.delta_count = delta_count;

%

lftfun = lft_finalize(lftfun);

%

end

Appendix C

141

C.2 Main Script

%% MAIN script – l’eseguibile che invoca tutte le altre

funzioni

clc, clear all, close all

% NB nell'errore l'errore "Index exceeds the number of array

elements (5).

% il numero (#) si riferisce ai valori normalizzati dei delta,

che si

% definiscono nel main

%

%% PARAMETERS TO FIND

f_vfa_true = 0.482254;

k_la_h2_true = 0.0004861;

k_la_true = 0.002305;

%

%parameter input

load_parameters;

%

%% (ACTIVE) unknown parameters range definition

f_vfa_lim=[0 1];

k_la_h2_lim=[0 10];

k_la_lim=[0 10];

%

%% (PASSIVE) data loading

% N.B. devono essere definiti vettori riga di condizioni

iniziali nonchè

% array di di variabili - dataseries di input, output e serie

temporale

load('data.mat');

load('initialcnd.mat');

%

% (ACTIVE) define the starting conditions)

% initialcnd is a table that contains the state initialization

% values

initial_cnd = [initialcnd.X_01_start,...

initialcnd.X_nd_start,...

initialcnd.X_1_start,...

initialcnd.X_2_start,...

initialcnd.X_3_start,...

initialcnd.S_1_start,...

initialcnd.S_2_start,...

initialcnd.S_4_start,...

initialcnd.S_5_start,...

initialcnd.IC_start,...

initialcnd.IN_start,...

initialcnd.S_hco3_start,...

initialcnd.S_nh3_start,...

Appendix C

142

initialcnd.S_cat_start,...

initialcnd.S_an_start,...

initialcnd.S_gas_h2_start,...

initialcnd.S_gas_ch4_start,...

initialcnd.S_gas_co2_start,...

]';

%

% (ACTIVE) define the input variables

in_u=[data.X_01_in,...

data.X_1_in,...

data.X_2_in,...

data.X_3_in,...

data.S_1_in,...

data.S_2_in,...

data.IC_in,...

data.IN_in,...

data.S_cat_in,...

data.S_an_in,...

];

%

% (ACTIVE) define the output variables

lft_y=[data.X_T,...

data.S_1,...

data.S_2,...

data.Alk,...

data.q_gas,...

data.S_gas_h2,...

data.S_gas_ch4,...

data.S_gas_co2,...

data.S_h

];

%//////////////////////////////////////////////

%% /////// OPTIONS /////////////////7///////////

%/////////////////////////////////////////////

LFTsolverOptions = lftSet('SensAlgorithm', 'ode15s',...

'SolutionTimeSpan', data.time,...

'SolutionInterpMethod', 'spline',...

'OversamplingMethod', 'spline',...

'RelTol', 1e-3,...

'AbsTol', 1e-6...

);

LFToptimOptions = lftOptSet_old('StartOptimSample', 1,...

'Display', 'iter',...

'StepTolerance', 1e-6,...

'MaxIter', 120, ...

'TolFun', 1e-9 ...

);

%

%% (PAssive) Input and output definition

Input = struct('Type', 'interpolated', 'Samples', in_u,

'Time', data.time);

Appendix C

143

InitialConditions = struct('StateInitialConditions',

initial_cnd);

%

%% LFT FUNCTION CALL

% must have the samee parameters of the lftfun function

lftfun = AMOCO_2_18_lft(R,...

T,...

V,...

V_gas,...

Q,...

p_gas,...

p_h2o,...

mu_1_max,...

mu_2_max,...

mu_3_max,...

k_d_1,...

k_d_2,...

k_d_3,...

K_s_1,...

K_s_2,...

K_s_4,...

Y_x_1,...

Y_x_2,...

Y_x_3,...

k_h1,...

k_nd,...

C_bac,...

C_x01,...

C_s1,...

C_s2,...

C_ch4,...

N_bac,...

N_x01,...

N_s1,...

k_ab_co2,...

k_a_co2,...

k_ab_IN,...

k_a_IN,...

k_w,...

k_h_h2,...

k_h_ch4,...

k_h_co2,...

f_vfa_lim,...

k_la_h2_lim,...

k_la_lim);

%

% (ACTIVE) same number of zeroes as the delta matrix dimension

lftfun.DeltaVal = diag([0 0 0 0 0 0 0 0 0 0 0 0 0 ]);

%

%% W matrix definition (PASSIVE)

M=max(lft_y);

Appendix C

144

m=min(lft_y);

for i=1:length(M)

nominal(i)=(M(i)-m(i));

end

lftfun.ISO.Nominal=nominal;

%

%% ESTIMATE (PASSIVE)

%

[DELTA_opt,fval,J,grad,H,CN,history] =

lftOptDelta(lftfun,Input,InitialConditions,LFTsolverOptions,lf

t_y,LFToptimOptions);

%

norm2abs(lftfun, lftfun.DeltaVal, DELTA_opt);

% displaying the real parameteric values

disp('f_h2_true = 0.482254');

disp('k_la_h2_true = 0.0004861');

disp('k_la_true = 0.002305');

%

H=H(:,:,end), CN=CN(end)

%% END

Appendix C

145

C.3 Load_Parameters script

This script is used to execute the parameter import

% parameters is a table that contains the parametric values

load('parameters.mat');

R = parameters.R ;

T = parameters.T ;

Q = parameters.Q ;

V_gas = parameters.V_gas ;

V = parameters.V ;

p_gas = 101325;

% BIOCHEICAL PARAMETERS

mu_1_max = parameters.mu_1_max ;

mu_2_max = parameters.mu_2_max ;

mu_3_max = parameters.mu_3_max ;

k_d_1 = parameters.k_d_1 ;

k_d_2 = parameters.k_d_2 ;

k_d_3 = parameters.k_d_3 ;

K_s_1 = parameters.K_s_1 ;

K_s_2 = parameters.K_s_2 ;

K_s_4 = parameters.K_s_4 ;

Y_x_1 = parameters.Y_x_1 ;

Y_x_2 = parameters.Y_x_2 ;

Y_x_3 = parameters.Y_x_3 ;

k_h1 = parameters.k_h1 ;

k_nd = parameters.k_nd ;

% NITROGEN, CARBON FRACTIONING

C_bac = parameters.C_bac ;

C_x01 = parameters.C_x01 ;

C_s1 = parameters.C_s1 ;

C_s2 = parameters.C_s2 ;

C_ch4 = parameters.C_ch4 ;

N_bac = parameters.N_bac ;

N_x01 = parameters.N_x01 ;

N_s1 = parameters.N_s1 ;

% CHEMICAL PARAMETERS

k_ab_co2= parameters.k_ab_co2 ;

k_ab_IN = parameters.k_ab_IN ;

k_a_co2 = parameters.k_a_co2 ;

k_a_IN = parameters.k_a_IN ;

k_w = parameters.k_w ;

p_h2o = parameters.p_h2o ;

k_h_co2 = parameters.k_h_co2 ;

k_h_ch4 = parameters.k_h_ch4 ;

k_h_h2 = parameters.k_h_h2 ;

% CALCULATED PARAMETERS

d = Q/V;

ro_1= k_h_h2*R*T*1000;

Appendix C

146

ro_2= k_h_ch4*R*T*1000;

ro_3= k_h_co2*R*T*1000;

ro_4= R*T/(p_gas-p_h2o)*V*1000;

ro_5= V/V_gas;

Appendix C

147

C.4 EquationsToMatrix script

The aim of this script is to aid the programmer in the creation of the LFT matrices starting

from the AMOCO openmodelica equations. The user must define the unknown paramteres

and guide the algorithm in the presence of denominator-bound non linearities.

clc

clear all

%% symbolic parameter loading

symbolic_parameters

%% DECLARE unknown parameters

unknown_parameters = [

f_vfa

k_la_h2

k_la

];

% DELTA matrix creation

syms DELTA_ [length(unknown_parameters),1]

%% DECLARE symbolic state variables

symbolic_state_variables

AMOCO_state_array = [ X_01

X_nd

X_1

X_2

X_3

S_1

S_2

S_4

S_5

IC

IN

S_hco3

S_nh3

S_cat

S_an

S_gas_h2

S_gas_ch4

S_gas_co2

];

% symbolic lft state variables creation

syms x_ [length(AMOCO_state_array),1]

%% DECLARE inputs

symbolic_input

AMOCO_input_array = [

X_01_in

Appendix C

148

X_1_in

X_2_in

X_3_in

S_1_in

S_2_in

IC_in

IN_in

S_cat_in

S_an_in

];

u = sym('u_%d',[length(AMOCO_input_array),1]);

%% other variables

syms r_ [1,13]

syms p_co2 p_h2 p_ch4

syms S_h

%% parameteri derivati

f_h2 = 1-f_vfa;

S_nh4 = IN - S_nh3;

S_3 = IC - S_hco3;

p_h2 = S_gas_h2*R*T;

p_ch4 = S_gas_ch4*R*T;

p_co2 = S_gas_co2*R*T;

%% rate equations

r_1 = k_h1 * X_01;

r_2 = mu_1_max / Y_x_1 * S_1 / (S_1 + K_s_1) * X_1;

r_3 = mu_2_max / Y_x_2 * S_2 / (S_2 + K_s_2) * X_2;

r_4 = mu_3_max / Y_x_3 * S_4 / (S_4 + K_s_4) * X_3;

r_5 = k_d_1 * X_1;

r_6 = k_d_2 * X_2;

r_7 = k_d_3 * X_3;

r_8 = k_la_h2 * (S_4 / 16 - k_h_h2 * p_h2) ;

r_9 = k_la * (S_5 / 64 - k_h_ch4 * p_ch4) ;

r_10 = k_la * (S_3 - k_h_co2 * p_co2) ;

r_12 = k_ab_co2 * (S_hco3 * (k_a_co2 + S_h) - k_a_co2 * IC);

r_13 = k_ab_IN * (S_nh3 * (k_a_IN + S_h) - k_a_IN * IN);

%% gas equation

q_gas = R * T / (p_gas - p_h2o) * V * 1000 * (r_8 + r_9 +

r_10);

%% state equations (copied verbatim)

X_01_dot = Q / V * (X_01_in - X_01) - r_1 + r_5 + r_6 + r_7;

X_nd_dot = Q / V * (-X_nd) + k_nd * r_1;

X_1_dot = Q / V * (X_1_in - X_1) + Y_x_1 * r_2 - r_5;

X_2_dot = Q / V * (X_2_in - X_2) + Y_x_2 * r_3 - r_6;

X_3_dot = Q / V * (X_3_in - X_3) + Y_x_3 * r_4 - r_7;

S_1_dot = Q / V * (S_1_in - S_1) + (1 - k_nd) * r_1 - r_2;

S_2_dot = Q / V * (S_2_in - S_2) + (1 - Y_x_1) * f_vfa *

r_2 - r_3;

Appendix C

149

S_4_dot = Q / V * ( - S_4) + (1 - Y_x_1) * f_h2 * r_2 - r_4

- r_8 * 16 ;

S_5_dot = Q / V * ( - S_5) + (1 - Y_x_2) * r_3 + (1 - Y_x_3)

* r_4 - r_9 * 64;

IC_dot = Q / V * (IC_in - IC)...

+ r_1 * (C_x01 - C_s1)...

+ r_2 * (C_s1 - (Y_x_1 * C_bac + (1 - Y_x_1) * C_s2 *

f_vfa))...

+ r_3 * (C_s2 - (Y_x_2 * C_bac + (1 - Y_x_2) * C_ch4))...

+ r_4 * (-(Y_x_3 * C_bac + (1 - Y_x_3) * C_ch4))...

+ (r_5 + r_6 + r_7) * (C_bac - C_x01)...

- r_10;

IN_dot = Q / V * (IN_in - IN)...

+ r_1 * (N_x01 - N_s1)...

+ r_2 * (N_s1 - Y_x_1 * N_bac)...

+ r_3 * (-Y_x_2 * N_bac)...

+ r_4 * (-Y_x_3 * N_bac)...

+ (r_5 + r_6 + r_7) * (N_bac - N_x01);

S_hco3_dot = -r_12;

S_nh3_dot = -r_13;

S_cat_dot = Q / V * (S_cat_in - S_cat) ;

S_an_dot = Q / V * (S_an_in - S_an);

S_gas_h2_dot = (-q_gas * S_gas_h2 / V_gas) + V / V_gas * r_8

;

S_gas_ch4_dot = (-q_gas * S_gas_ch4 / V_gas) + V / V_gas *

r_9;

S_gas_co2_dot = (-q_gas * S_gas_co2 / V_gas) + V / V_gas *

r_10;

%% equation array creation

eqns = [

X_01_dot

X_nd_dot

X_1_dot

X_2_dot

X_3_dot

S_1_dot

S_2_dot

S_4_dot

S_5_dot

IC_dot

IN_dot

S_hco3_dot

S_nh3_dot

S_cat_dot

S_an_dot

S_gas_h2_dot

S_gas_ch4_dot

S_gas_co2_dot

];

%% gas transfer equation substitution

Appendix C

150

G_array =[

(S_4/16 - k_h_h2*p_h2)

(S_5/64 - k_h_ch4*p_ch4)

(S_3 - k_h_co2*p_co2)

];

syms G_ [length(G_array),1]

eqns=subs(eqns,G_array,G);

%% MONOD substitution (monod)

%syms x_ [1 length(x)]

MONOD_array = [

S_1 / (S_1 + K_s_1) * X_1

S_2 / (S_2 + K_s_2) * X_2

S_4 / (S_4 + K_s_4) * X_3

];

syms MONOD_ [length(MONOD_array),1]

eqns =subs(eqns, MONOD_array,MONOD);

%% x u DELTA substitution

eqns_subs = subs(eqns,AMOCO_state_array,x);

eqns_subs = subs(eqns_subs, AMOCO_input_array, u);

eqns_subs = subs(eqns_subs, unknown_parameters, DELTA);

%% expansion

eqns_subs = expand(eqns_subs);

%% XI substitution

XI_array =[

MONOD_1

MONOD_2

MONOD_3

S_h*x_12

S_h*x_13

G_1*x_16

G_1*x_17

G_1*x_18

G_2*x_16

G_2*x_17

G_2*x_18

G_3*x_16

G_3*x_17

G_3*x_18

];

syms XI_ [length(XI_array),1]

eqns_final=subs(eqns_subs,XI_array,XI);

%% sostituzione delle w

%la matrice delle w completa è

w_array = [

DELTA_1*XI_1

DELTA_2*G_1

Appendix C

151

DELTA_2*XI_6

DELTA_2*XI_7

DELTA_2*XI_8

DELTA_3*G_2

DELTA_3*G_3

DELTA_3*XI_9

DELTA_3*XI_10

DELTA_3*XI_11

DELTA_3*XI_12

DELTA_3*XI_13

DELTA_3*XI_14

];

syms w_ [length(w_array),1]

eqns_final_w =subs(eqns_final,w_array,w);

%% matrice delle z:

z_array = subs(w_array, DELTA, [1;1;1;]);

z_array = subs(z_array, G, G_array);

z_array = subs(z_array,AMOCO_state_array,x);

%% matrice delle om

THETA = S_cat + S_nh4 - S_hco3 - S_an - S_2 / 64;

THETA_subs = subs(THETA,AMOCO_state_array,x);

S_h_subs = 1/2*(-THETA + sqrt(THETA^2 + 4*k_w));

% sostituzione delle variabili G, MONOD, S_h entro la matrice

delle XI

om_array = subs(XI_array, G, G_array);

om_array = subs(om_array, MONOD, MONOD_array);

om_array = subs(om_array, S_h, S_h_subs);

om_array = subs(om_array, AMOCO_state_array, x);

%% gas subs

syms k_1

eqns_gas = [

(k_la_h2*(G_1*V - G_1*S_gas_h2*k_1) -

k_la*(G_2*S_gas_h2*k_1 - G_3*S_gas_h2*k_1))/V_gas

(k_la*(G_2*V - G_2*S_gas_ch4*k_1 - G_3*S_gas_ch4*k_1) -

G_1*S_gas_ch4*k_1*k_la_h2)/V_gas

(k_la*(G_3*V - G_2*S_gas_co2*k_1 - G_3*S_gas_co2*k_1) -

G_1*S_gas_co2*k_1*k_la_h2)/V_gas

];

%%

eqns_gas_subs = subs(eqns_gas,unknown_parameters,DELTA);

eqns_gas_subs = subs(eqns_gas_subs, AMOCO_state_array, x);

%%

eqns_gas_subs_expand = expand(eqns_gas_subs);

eqns_gas_subs_expand = collect(eqns_gas_subs_expand,V_gas);

%%

array_gas_2 = [

(DELTA_2*G_1*V - DELTA_2*G_1*k_1*x_16 -

DELTA_3*(k_1*x_16*(G_2 + G_3)))/V_gas

Appendix C

152

(DELTA_3*G_2*V - DELTA_2*G_1*k_1*x_17 -

DELTA_3*(k_1*x_17*(G_2 + G_3)))/V_gas

(DELTA_3*G_3*V - DELTA_2*G_1*k_1*x_18 -

DELTA_3*(k_1*x_18*(G_2 + G_3)))/V_gas

];

%% XI substitution

XI_array =[

MONOD_1

MONOD_2

MONOD_3

S_h*x_11

S_h*x_13

G_1*x_16

G_1*x_17

G_1*x_18

x_16*(G_2 + G_3)

x_17*(G_2 + G_3)

x_18*(G_2 + G_3)

];

syms XI_ [length(XI_array),1]

eqns_gas_XI=subs(array_gas_2,XI_array,XI);

%%

w_gas = [

DELTA_1*XI_1

DELTA_2*G_1

DELTA_2*XI_6

DELTA_2*XI_7

DELTA_2*XI_8

DELTA_3*G_2

DELTA_3*G_3

DELTA_3*XI_9

DELTA_3*XI_10

DELTA_3*XI_11

];

%%

syms w_ [length(w_gas),1]

eqns_gas_w = subs(eqns_gas_XI,w_gas,w);

%% matrice delle z:

z_gas = subs(w_gas, DELTA, [1;1;1;]);

z_gas = subs(z_gas, G, G_array);

z_gas = subs(z_gas,AMOCO_state_array,x);