Energy balance: Difference between revisions

From Opasnet
Jump to navigation Jump to search
(OpasnetUtils::combine)
 
(7 intermediate revisions by 2 users not shown)
Line 1: Line 1:
<noinclude>
[[op_fi:Energiatase]]
[[op_fi:Energiatase]]
[[Category:Energy]]
[[Category:Energy]]
Line 6: Line 7:
[[Category:Code under inspection]]
[[Category:Code under inspection]]
{{method|moderator=Jouni|stub=Yes}}
{{method|moderator=Jouni|stub=Yes}}
</noinclude>


{{summary box |
question = What is energy balance and how is it modelled? |
answer =
Summing up the amount of energy produced and subtracting the amount of energy consumed within a time period gives the energy balance. Since the electricity grid and district heat network lack significant storage mechanics, the balance has to be virtually zero over short periods. When considering the balance of a particular area (e.g. Helsinki), we can make the assumption that electricity can be imported and exported in international markets. The energy in the district heat network, however, has to produced locally. This sets up the non-trivial problem of optimising production so that there are no significant deficits as well as minimising losses and maximising profits. This problem is solved (to some extent) by market forces in the real world. Our most recent energy balance model uses linear programming tools to solve an optimum for the activity of a given set of production units in simulated instances created by the main model. The main model is responsible for the decision making aspects, while the energy balance optimisation only functions as an approximation of real world market mechanics.
The linear programming problem is set up as follows.
For each production unit: let x<sub>i</sub> be activity of the plant. Lets also have variables y<sub>j</sub> for deficits and excesses for each type of energy produced.
The objective function is the function we are optimising. Each production unit has a unit profit per activity denoted by a<sub>i</sub> which is determined by the amount of different input commodities (e.g. coal) per amount of different output commodities (i.e. electricity and heat) and their market prices. Also, lets say we want to make sure that district heat demand is always met when possible and have a large penalty factor for each unit of heat demand not met (1M€ in the model). In addition, it must be noted that excess district heat becomes wasted so it counts as loss. Let these deficit and excess related losses be denoted by b<sub>j</sub>. The whole objective function then becomes: sum(x<sub>i</sub>a<sub>i</sub>) + sum(y<sub>j</sub>b<sub>j</sub>).
The values of variables are constrained by equalities and inequalities: the sum of production of a commodity is equal to its demand minus deficit plus excess, activity is constrained by the maximum capacity and all variables are non-negative by definition. 
This can be efficiently solved by computers for each given instance. Production wind-up and wind-down is ignored, since time continuity is not considered. As a consequence fuel limits (e.g. diminishing hydropower capacity) are not modelled completely either.
}}
==Question==
==Question==


How to calculate energy balances?
What is energy balance and how is it modelled?  


==Answer==
==Answer==


Press button to run the model. Get a free user account and edit the table below to change the model inputs. You may also copy the model to your own page for your own purposes. For examples, see e.g. [[Energy balance in Kuopio]] or [[Energy balance in Suzhou]].
Summing up the amount of energy produced and subtracting the amount of energy consumed within a time period gives the energy balance. Since the electricity grid and district heat network lack significant storage mechanics, the balance has to be virtually zero over short periods. When considering the balance of a particular area (e.g. Helsinki), we can make the assumption that electricity can be imported and exported in international markets. The energy in the district heat network, however, has to be produced locally. This sets up the non-trivial problem of optimising production so that there are no significant deficits as well as minimising losses and maximising profits. This problem is solved (to some extent) by market forces in the real world.  


<rcode name="answer" graphics="1">
In Opasnet, there are two different ways to calculate energy balance. Our most recent energy balance model uses linear programming tools to solve an optimum for the activity of a given set of production units in simulated instances created by the main model. The main model is responsible for the decision making aspects, while the energy balance optimisation only functions as an approximation of real world market mechanics. This version was used in [[Helsinki energy decision 2015]].
library(OpasnetUtils)
library(ggplot2)


objects.latest("Op_en5141", code_name = "initiate")
The previous version was based on setting up a set of linear equations describing the inputs, outputs, and shares of different energy and plant processes. This approach is less flexible, because it does not use an optimising function and everything must be described as linear (or piecewise linear). However, this approach was successfully used in [[Energy balance in Kuopio]] and [[Energy balance in Suzhou]].


N <- 10
<gallery heights="500" widths="500">
File:District heat flow in Helsinki.png|Example of energy balance model: District heat flow in Helsinki. The scenarios are from the assessment [[Helsinki energy decision 2015]].
File:Flow of electric power in Helsinki.png|Example of energy balance model: Electric power heat flow in Helsinki. The scenarios are from the assessment [[Helsinki energy decision 2015]].
File:Incomes and costs of energy production in Helsinki.png|Example of energy balance model: Incomes and costs of energy production in Helsinki. The scenarios are from the assessment [[Helsinki energy decision 2015]].
</gallery>


balance <- Ovariable("balance", ddata = "Op_en5141.equations")
This code is an example how the energy balance model is used in a city case. The data comes from [[Helsinki energy decision 2015]].


balance@data$Equation <- as.character(levels(balance@data$Equation)[balance@data$Equation]) # Change factor into character.
<rcode name="answer" graphics="1">
balance@data$Policy[balance@data$Policy == ""] <- NA # Prepare indices for fillna. This should be part of fillna.
library(OpasnetUtils)
balance@data$CHPcapacity[balance@data$CHPcapacity == ""] <- NA
library(ggplot2)
balance@data <- fillna(balance@data, marginals = 1:2) # Fill empty slots in indices.


nonlinearity <- Ovariable("nonlinearity", ddata = "Op_en5141", subset = "Nonlinearity parameters")
openv.setN(0) # use medians instead of whole sampled distributions


directinput <- Ovariable("directinput", ddata = "Op_en5141", subset = "No modelled upstream variables")
# Download case-specific data, in this case from Helsinki.
objects.latest("Op_en7237", code_name = "intermediates") # [[Helsinki energy decision 2015]]


# Download energy balance model and its parts:
#  EnergyConsumerDemandTotal
#  EnergyFlowHeat
#  EnergyFlowOther
#  EnergyNetworkDemand
#  EnergyNetworkOptim
#  fuelUse
#  EnergyNetworkCost
objects.latest("Op_en5141", code_name = "EnergyNetworkOptim") # [[Energy balance]]


eb <- EvalOutput(energy.balance)
oprint(summary(EvalOutput(EnergyNetworkOptim)))
 
eb@marginal[colnames(eb@output) == "energy.balanceVar"] <- TRUE # This has to be done manually because CheckMarginals does not notice this as a marginal.
eb@marginal[colnames(eb@output) %in% nonlinearity@output$critIndex] <- FALSE # Nonlinear indices are demarginalised because only one of the two equations apply.
 
oprint(summary(eb))
 
ggplot(eb@output, aes_string(x = "energy.balanceVar", y = "energy.balanceResult", fill = "Policy")) +
geom_boxplot() +
theme_grey(base_size = 24) +
theme(axis.text.x = element_text(angle = 90, hjust = 1))


</rcode>
</rcode>
Line 62: Line 54:
==Rationale==
==Rationale==


===Input===
=== Energy balance with linear programming ===


; Example table for making matrices from text format equations. CHPcapacity describes which of the piecewise linear equations should be used. Policy is a decision option that alters the outcome. Dummy is only for compatibility but it is not used.
The linear programming problem is set up as follows.


<t2b name="Equations" index="CHPcapacity,Policy,Equation" obs="Dummy" desc="Description" unit="GWh /a">
For each production unit: let x<sub>i</sub> be activity of the plant. Lets also have variables y<sub>j</sub> for deficits and excesses for each type of energy produced.
|Biofuel|CHP renewable = CHP peat|1|Biofuel policy contains half biofuels, half peat
|BAU|CHP renewable = 89.24|1|
||CHP peat + CHP renewable + CHP oil = CHP heat + CHP electricity|1|
||CHP peat = 90-98*CHP oil|1|
||CHP electricity = 0.689*CHP heat|1|
CHP<1000||H heat = 0.08*CHP heat|1|Small heat plants reflect the total heat need
CHP>1000||CHP heat + CHP electricity = 1000|1|But production capacity of CHP may be overwhelmed, decoupling CHP heat and H heat.
||H biogas + H oil = H heat|1|
||H oil = 18.973*H biogas|1|
||Bought electricity + CHP electricity = Cons electricity|1|
||CHP heat + H heat = Cons heat|1|
||Cons electricity = 900-1100|1|
||Cons heat = 900-1000|1|
</t2b>


; Example table to describe the details about nonlinear equations.
The objective function is the function we are optimising. Each production unit has a unit profit per activity denoted by a<sub>i</sub> which is determined by the amount of different input commodities (e.g. coal) per amount of different output commodities (i.e. electricity and heat) and their market prices. Also, lets say we want to make sure that district heat demand is always met when possible and have a large penalty factor for each unit of heat demand not met (1 M€ in the model). In addition, it must be noted that excess district heat becomes wasted so it counts as loss. Let these deficit and excess related losses be denoted by b<sub>j</sub>. The whole objective function then becomes: sum(x<sub>i</sub>a<sub>i</sub>) + sum(y<sub>j</sub>b<sub>j</sub>).  


<t2b name="Nonlinearity parameters" index="critVar,critIndex,rescol,critLocLow,critLocHigh" obs="critValue" unit="GWh /a">
The values of variables are constrained by equalities and inequalities: the sum of production of a commodity is equal to its demand minus deficit plus excess, activity is constrained by the maximum capacity and all variables are non-negative by definition.
Cons heat|CHPcapacity|Result|CHP<1000|CHP>1000|1080
This can be efficiently solved by computers for each given instance. Production wind-up and wind-down is ignored, since time continuity is not considered. As a consequence fuel limits (e.g. diminishing hydropower capacity) are not modelled completely either.
</t2b>
 
 
; This table is fetched if there are no nonlinearities. Therefore, there is no need to copy it to the case study page.
 
<t2b name="No nonlinearities" index="critVar,critIndex,rescol,critLocLow,critLocHigh" obs="critValue" unit="GWh /a">
|||||
</t2b>
 
; This table is fetched if there are no modelled upstream variables that would affect the equations.
 
<t2b name="No modelled upstream variables" index="energybalanceVars" obs="Result" unit="-">
|
</t2b>
 
===Output===
 
Output is an ovariable with column energy.balanceVar for the solved variables, a result column, and all other indices in the balance object used as input. This method has been used in several cases:
* [[Energy balance in Kuopio]]
* [[Energy balance in Basel]]
* [[Energy balance in Stuttgart]]
* [[Energy balance in Suzhou]]
 
===Data===
 
* Energy balances are described as input = output on a coarse level (called classes) where the structure is the same or similar to the OECD energy balance tables. If possible, this is described on the Energy balance method level and it is shared by all cities.
* On more detailed (variable level in the matrix), the fraction of each variable of the total class are described separately. Fractions are city specific and they are described on city level in a separate table.
* Based on the fraction table, detailed equations with variables are created. The format will be fraction * class total = variable.
* The last fraction has zero degrees of freedom when the class total is given. However, it must have a variable and thus a row in the fraction table. The result for that variable is an empty cell (which results in NA).
* Unlike in the previous version, all variables are given either as values or equations, and the user interface is not used for BAU. In contrast, user interface or decision table may be used to derive values for alternative scenarios.
* To make this work, the city-specific fraction data must be defined as ovariable (so that it can be changed with a decision table), and also the energy balance method must be described asa ovariable. How are we going to make the two interplay, as we may want to have several cities?
** Define one city ovariable and evaluate energy balance with that. The ovariable has a generic name. Then, define a new city ovariable with the same name and re-evaluate the energy balance ovariable; this must be done so that the two cities are appended rather than replaced.
** city ovariables are appended first into a large fraction table, and then that is used to create the large energy balance matrix. {{defend|# |This is clearly better.|--[[User:Jouni|Jouni]] 17:09, 21 February 2013 (EET)}}
* The city-specific ovariable may have Iter and other indices. A separate matrix is created and solved for each unique combination of indices. This makes it possible to have a very flexible approach.
* We should check if the energy balance matrix (see Matti's Excel) has city-specific equations. If possible, energy transformations are described as generic equations on the energy balance method.
* Structure of OECD Energy balance tables (data):
** Fuel (given as observation columns in OECD table)
** Activity (row in OECD table)
** Description
* Structure of the generic process table
** Equation,
** Col,
** Result,
** Description? {{attack|# |This does not join up in a coherent way.|--[[User:Jouni|Jouni]] 17:09, 21 February 2013 (EET)}}
* Columns for fraction table
** Class
** Item
** Result (fraction)
** Indices as needed


''''Additional thoughts
:''Ovariables like EnergyNetworkOptim below are used in [[Helsinki energy decision 2015]]. [[Prices of fuels in heat production]] are used as direct inputs in the optimising.
 
Energy balance has these parts
* General table that describes balances and conversions from ole fuel to another.
** Items and sums Are defined at generic level.
* Data table that contains amounts from OECD table
* Fraction table thAt describes how detailSs Are derived from sums.
* if there Are additions to sta dard items, they descrided at city level in an addition table. It must have the struxture of the matrix table in long format.
* can nonlinearities ne handled simply by indices separating two linear parts of nonlinear equations?
 
====Old description====
 
{{attack|# |What are other energy sources,i wish it should be specific,how about air transport,does it included in the energy balance calculation?|--[[User:Sam0911|Sam0911]] 16:29, 10 February 2013 (EET)}}
 
Use a table where different fuel types are columns and different stocks, energy production processes, or consumption types are rows called ''Energy accounts'' (Ene.account for short). See also an example [[File:Energy supply in Europe.xls]]. All Ene.accounts and fuel types used are listed on [[Energy consumption classes]].
 
{| {{prettytable}}
|+ '''A part of a full energy balance sheet (ktoe/a).
|----
!
! Coal and peat
! Crude oil
! Petrochemical products
! Electricity
! Heat
|----
!colspan="6"|Production and supply of primary energy
|----
| Production || 129|| || || ||
|----
| Import || || 28|| 63|| ||
|----
| Export || || || || || 
|----
!colspan="6"| Conversion of primary energy to use
|----
| Transfers || || || || ||
|----
| Electricity plants || || || || ||
|----
| CHP plants || -129|| -28|| || 61|| 96
|----
!colspan="6"|Final energy consumption
|----
| Industry || || || || ||
|----
| Road transport || || || 54|| ||
|----
| Heating || || || || 10|| 96
|----
| Other energy use || || || 9|| 51||
|}
 
Energy balance can be considered as double-entry bookkeeping. The production and conversion are typically credit (i.e., where the energy comes from), while consumption is debit (i.e., where the energy is used). In a typical case, both credit and debit are marked positive on the energy balance sheet. However, there are activities that convert one fuel to another, so that the energy is moved from one column to another. In this case, credit (the source) is marked negative and debit (the target) is marked positive. This is because the production and conversion rows together should reflect how much energy is actually available for final consumption. In other words, the sum of production and conversion rows should equal the sum of the consumption rows in every column.
 
Because production + conversion = consumption, also production = consumption - conversion. These equations are used to derive the supply that is needed to fulfil the demand.
 
=== Fuel use and fuel shares in generic processes ===
 
There is an alternative to energy for calculating fuel use. It is based on the idea that for each heating type, there is a constant share of fuels used. For some heating types, this is generic and is shown on this page. For some others, the constant is case-specific and is determined on a case-specific page.
 
The table below contains connections of heating types and fuel usage in generic situations. There may be case-specific differences, which must be handled separately.
 
<t2b name='Fuel use in different heating types' index='Heating,Burner,Fuel' obs='Fraction' desc='Description' unit='-'>
Wood|Household|Wood|1|
Oil|Household|Light oil|1|
Gas|Household|Gas|1|
Heating oil|Household|Light oil|1|
Coal|Household|Coal|1|
Other sources|Household|Other sources|1|
No energy source|Household|Other sources|1|
Geothermal|Grid|Electricity|0.3|Geothermal does not sum up to 1 because more heat is produced than electricity consumed.
Centrifuge, hydro-extractor|Grid|Electricity|0.3|Not quite clear what this is but presumably a heat pump.
Solar heater/ collector|Grid|Electricity|0.1|Use only; life-cycle impacts omitted.
Electricity|Grid|Electricity|1|
District|Undefined|Heat|1|
</t2b>
{{comment|# |District removed from fuelShares/Heating, but it is not clear why it was there in the first place. For completion? But for what?|--[[User:Jouni|Jouni]] ([[User talk:Jouni|talk]]) 13:41, 23 July 2015 (UTC)}}
:{{attack|# |Still necessary in model. It completes the description of flows related to different heating methods. In a descriptor of flows, "fuel" is a slight misnomer and should be reclassified as a "commodity" or similar throughout the model.|--[[User:Teemu R|Teemu R]] ([[User talk:Teemu R|talk]]) 12:23, 24 July 2015 (UTC)}}
 
<rcode name='fuelSharesgeneric' embed=1 label='Initiate fuelSharesgeneric (only for developers)'>
 
# This is code Op_en5141/fuelSharesgeneric (only generic) on page [[Energy balance]].
library(OpasnetUtils)
 
fuelSharesgeneric <- Ovariable("fuelSharesgeneric", ddata = "Op_en5141", subset = "Fuel use in different heating types") # [[Energy balance]]
colnames(fuelSharesgeneric@data) <- gsub("[ \\.]", "_", colnames(fuelSharesgeneric@data))
#fuelSharesgeneric@data <- merge(fuelSharesgeneric@data, data.frame(Time = 1900:2080))
 
objects.store(fuelSharesgeneric)
cat("Object fuelSharesgeneric initiated!\n")
 
</rcode>
 
 
<rcode name="fuelUse" label="Initiate fuelUse (developers only)" embed=1>
## This is code Op_en5141/fuelUse on page [[Energy balance]].
library(OpasnetUtils)
 
fuelUse <- Ovariable("fuelUse",
dependencies = data.frame(
Name = c("energyUse", "fuelShares")
),
formula = function(...) {
out <- energyUse * fuelShares
return(out)
}
)
 
objects.store(fuelUse)
cat("Ovariable fuelUse stored.\n")
</rcode>
 
=== Calculations ===
 
:''Ovariable energyBalance below is used in [[Helsinki energy decision 2015]]. [[Prices of fuels in heat production]] are used as direct inputs in the optimising.
 
<rcode name="energyBalancetest" label="Initiate energyBalance (developers only)" embed=1>
## This code is  Op_en5141/energyBalance [[Energy balance]].
library(OpasnetUtils)
 
energyBalance <- Ovariable("energyBalance",
dependencies = data.frame(
Name = c(
"energy", # Energy supply and demand in the system
"nondynsupply", # Energy supply that is not optimised.
"requiredName", # Fuel name for the energy type that is balanced for inputs and outputs.
"fuelShares", # Shares of fuels for different heating types.
"energyProcess", # Matrix showing inputs and outputs of the process of each plant.
"fuelPrice", # Prices of fuels (by Time)
"timelylimit", # Maximum energy supply in the given conditions at a particular time point
"plantParameters" # Capacities and costs of running each plant
),
Ident = c(
NA,
NA,
NA,
NA,
NA,
"Op_en4151/fuelPrice",
NA,
NA
)
),
formula = function(...) {
# First, define the optimising function
optfunction <- function(
x, # Variables to optimise (activities of power plants).
demand, # energy demand that must be met by supply (scalar)
unitprice, # vector of prices per activity unit (by plant).
unitreq # vector of required energy produced per activity uniit (by plant).
) {
 
out <- sum(x * unitprice) + (demand - sum(x * unitreq))^2 * 1E6
return(out)
}
 
# Simplify the energy demand and supply data and add index Fuel.
ene <- unkeep(energy, prevresults = TRUE, sources = TRUE)
ene <- oapply(ene, cols = c("Efficiency", "Renovation", "Building"), FUN = sum) # I: X, Heating, Consumamble, Fuel
 
# Take the heating and look for fuel shares.
eneheating <- unkeep(ene, cols = "Fuel")
eneheating@output <- eneheating@output[eneheating@output$Heating != "Not heating" , ]
eneheating <- oapply(eneheating * fuelShares, cols = "Heating", FUN = sum) # I: X, Consumable + Burner, Fuel
 
# Then take the non-heating, which already has fuel shares.
enenonheating <- ene * Ovariable(output = data.frame(Burner = "Unidentified", Result = 1), marginal = c(TRUE, FALSE))
enenonheating@output <- enenonheating@output[enenonheating@output$Heating == "Not heating" , ]
enenonheating <- unkeep(enenonheating, cols = "Heating") # I: X, Consumable, Fuel + Burner
 
ene <- combine(eneheating, enenonheating) # I: X, Consumable, Fuel, Burner
 
if(!is.na(requiredName)) { # Energy balance is only done if requiredName is given
 
# Non-dynamic energy supply and demand, which is not optimised.
# It is multiplied by -1 because energy production consumes fuels.
nondynen <- ene
nondynen@output <- nondynen@output[nondynen@output$Fuel != requiredName , ]
nondynen <- nondynen * Ovariable( # I: X, Consumable, Fuel, Burner, Plant
output = data.frame(Plant = "Unidentified", Result = -1),
marginal = c(TRUE, FALSE)
)
 
# Dynamic energy supply and demand, which IS optmised.
dynen <- ene
dynen@output <- dynen@output[dynen@output$Fuel == requiredName , ]
dynen <- oapply(dynen, cols = c("Fuel", "Burner", "Consumable"), FUN = sum) # I: X
 
# Operational costs per unit of activity. Fuel price is excluded and calculated elsewhere.
operationCost <- plantParameters
operationCost@output <- operationCost@output[operationCost@output$Parameter == "Operation cost" , ]
operationCost <- unkeep(operationCost, cols = "Parameter") # I: Y, Plant
 
# Costs per unit of activity by plant. Fuel price included. Also revenues included.
unitprice <- oapply(energyProcess * fuelPrice * -1, cols = "Fuel", FUN = sum) + operationCost
colnames(unitprice@output)[colnames(unitprice@output) == "plantParametersResult"] <- "unitpriceResult"
unitprice@output$unitpriceResult <- result(unitprice) # Save result into another column. I: Z, Plant, Burner - Fuel
# Required energy produced per unit activity in each plant. Only take plants that affect required energy.
unitreq <- energyProcess
unitreq@output <- unitreq@output[unitreq@output$Fuel == requiredName & result(unitreq) != 0 , ]
unitreq <- unkeep(unitreq, cols = "Fuel") # I: Plant, Burner - Fuel
# Temporary limit to energy supply due to time-related particular constraints. Exclude unused plants.
timlim <- plantParameters
timlim@output <- timlim@output[timlim@output$Parameter == "Max" , ]
timlim <- (
unkeep((timelylimit >= timlim) * timlim, prevresults = TRUE, sources = TRUE) +
(timelylimit < timlim) * timelylimit
)
timlim@output <- timlim@output[result(timlim) != 0 , ] # Exclude if max == 0
colnames(timlim@output)[colnames(timlim@output) == "plantParametersResult"] <- "timelylimitResult"
timlim@output$timelylimitResult <- result(timlim) # Save result into another column
timlim <- unkeep(timlim, cols = "Parameter") # I: Y, Plant
# Take the minimum from the plantParameters.
ppar <- plantParameters
ppar@output <- ppar@output[ppar@output$Parameter == "Min" , ]
ppar <- unkeep(ppar, cols = "Parameter", sources = TRUE, prevresults = TRUE)
 
# Put all inputs into one ovariable and sort by Plant. Note: important outputs in columns
# unitpriceResult, energyProcessResult, plantParametersResult, timelylimitResult and Result (dynen).
temp <- unitprice * unitreq * timlim * ppar * 0 + dynen # I: X + Y + Z, Plant, Burner
 
indices <- unique(temp@output[temp@marginal & !colnames(temp@output) %in% c("Plant", "Burner", "Parameter")])
 
out <- data.frame()
for(i in 1:nrow(indices)) {
 
temp2 <- merge(indices[i , , drop = FALSE], temp@output)
temp2 <- temp2[order(temp2$Plant) , ]
temp3 <- optim(
par = temp2$plantParametersResult,
fn = optfunction,
lower = temp2$plantParametersResult,
upper = temp2$timelylimitResult,
unitprice = temp2$unitpriceResult,
demand = temp2$Result[1], # from dynen
unitreq = temp2$energyProcessResult,
method = "L-BFGS-B"
)
 
out <- rbind(out, data.frame(
indices[i , , drop = TRUE],
Plant = temp2$Plant,
Result = temp3$par
))
}
out <- Ovariable(output = out, marginal = c(rep(TRUE, ncol(out)-1), FALSE))
out <- out * energyProcess # Optimised fuel usage
out <- combine(out, nondynen) # Combine optimised and non-optimised energy demand
} else {
out <- ene # If optimisation is not done, take the simple version updated with fuelShares.
}
out <- unkeep(out, source = TRUE)
marginals <- colnames(out@output)[out@marginal]
out <- combine(out, nondynsupply) # Combine also non-dynamic energy supply
# Fill in other indices (typically decisions) that are missing in nondynsupply
marginals <- setdiff(marginals, colnames(nondynsupply@output)[nondynsupply@marginal])
out@output <- fillna(out@output, marginals = match(marginals, colnames(out@output)))
 
return(out)
}
)
 
objects.store(energyBalance)
cat("Ovariable energyBalance stored.\n")
</rcode>
 
==== Newer version ====


<rcode name="EnergyNetworkOptim" label="Initiate EnergyNetworkOptim(developers only)" embed=1>
<rcode name="EnergyNetworkOptim" label="Initiate EnergyNetworkOptim(developers only)" embed=1>
Line 454: Line 108:
EnergyFlowOther <- merge(
EnergyFlowOther <- merge(
EnergyConsumerDemandTotal[!EnergyConsumerDemandTotal$Consumable %in% c("Heating", "Hot water"),],  
EnergyConsumerDemandTotal[!EnergyConsumerDemandTotal$Consumable %in% c("Heating", "Hot water"),],  
unique(combine(temperene, nontemperene)[,c("Consumable", "Fuel")])
unique(OpasnetUtils::combine(temperene, nontemperene)[,c("Consumable", "Fuel")])
)
)
EnergyFlowOther@name <- EnergyConsumerDemandTotal@name
EnergyFlowOther@name <- EnergyConsumerDemandTotal@name
Line 470: Line 124:
),
),
formula = function(...) {
formula = function(...) {
demand <- combine(
demand <- OpasnetUtils::combine(
EnergyFlowHeat[EnergyFlowHeat$Fuel %in% c("Heat", "Electricity"), colnames(EnergyFlowHeat@output) != "Burner"],
EnergyFlowHeat[EnergyFlowHeat$Fuel %in% c("Heat", "Electricity"), colnames(EnergyFlowHeat@output) != "Burner"],
EnergyFlowOther[EnergyFlowOther$Fuel %in% c("Cooling", "Electricity"),]
EnergyFlowOther[EnergyFlowOther$Fuel %in% c("Cooling", "Electricity"),]
Line 601: Line 255:
optf <- function(varsi) {
optf <- function(varsi) {
# Acquire relevant sections of variables with respect to index iteration
# Acquire relevant sections of variables with respect to index iteration
Line 736: Line 390:
return(outi)
return(outi)
}
}
#total_time <- Sys.time() - total_time
#aggregate(out$Result[out$Process_variable_name == "Profit"], out[out$Process_variable_name == "Profit",colnames(out) %in% c("Time", "Temperature") | grepl("Policy", colnames(out))], mean)
#aggregate(out$Result[out$Process_variable_name == "Profit"], out[out$Process_variable_name == "Profit",colnames(out) %in% c("Time", "Temperature") | grepl("Policy", colnames(out))], mean)
#aggregate(out$Result[out$Process_variable_name == "HeatDeficit"], out[out$Process_variable_name == "HeatDeficit", colnames(out) %in% c("Time", "Temperature") | grepl("Policy", colnames(out))], mean)
#aggregate(out$Result[out$Process_variable_name == "HeatDeficit"], out[out$Process_variable_name == "HeatDeficit", colnames(out) %in% c("Time", "Temperature") | grepl("Policy", colnames(out))], mean)
#aggregate(out$Result[out$Process_variable_name == "Hanasaari"], out[out$Process_variable_name == "Hanasaari", colnames(out) %in% c("Time", "Temperature") | grepl("Policy", colnames(out))], mean)
#aggregate(out$Result[out$Process_variable_name == "Hanasaari"], out[out$Process_variable_name == "Hanasaari", colnames(out) %in% c("Time", "Temperature") | grepl("Policy", colnames(out))], mean)
out <- ddply(vars@output, colnames(vars@output)[vars@marginal & !colnames(vars@output) %in% exclude], optf)
# Error handling
optfsecure <- function(varsi) {
ret <- tryCatch(
optf(varsi),
error = function(e) return(NULL)
)
if (!is.null(ret)) {
return(ret)
} else {
warning(paste("EnergyNetworkOptim failed optimising a permutation with error:", geterrmessage()))
return(data.frame())
}
}
out <- ddply(vars@output, colnames(vars@output)[vars@marginal & !colnames(vars@output) %in% exclude], optfsecure)
out <- Ovariable(output = out, marginal = !grepl("Result$", colnames(out)))
out <- Ovariable(output = out, marginal = !grepl("Result$", colnames(out)))
Line 802: Line 469:
heating$Heating <- NULL
heating$Heating <- NULL
EnergyFlow <- combine(0 - EnergyFlow, real, heating)
EnergyFlow <- OpasnetUtils::combine(0 - EnergyFlow, real, heating)
EnergyFlow <- unkeep(EnergyFlow, sources = TRUE, prevresults = TRUE)
EnergyFlow <- unkeep(EnergyFlow, sources = TRUE, prevresults = TRUE)
#EnergyFlowTest <- oapply(EnergyFlow, c("Time", "Temperature", "Fuel"), sum)
#EnergyFlowTest <- oapply(EnergyFlow, c("Time", "Temperature", "Fuel"), sum)
Line 847: Line 514:


operation <- EnergyNetworkOptim[EnergyNetworkOptim@output$Process_variable_name == "Operation cost" , ]
operation <- EnergyNetworkOptim[EnergyNetworkOptim@output$Process_variable_name == "Operation cost" , ]
operation <- operation * temperdays * 5 * 1E-6 # For 5-year periods, € -> M€
operation <- operation * temperdays * 10 * 1E-6 # For 10-year periods, € -> M€
operation <- oapply(operation, cols = c("Temperature"), FUN = sum)
operation <- oapply(operation, cols = c("Temperature"), FUN = sum)
operation <- unkeep(operation, cols = c("Process_variable_name", "Process_variable_type"), sources = TRUE, prevresults = TRUE)
operation <- unkeep(operation, cols = c("Process_variable_name", "Process_variable_type"), sources = TRUE, prevresults = TRUE)
operation <- operation * Ovariable(output = data.frame(Plant = "Operation", Result = 1), marginal = c(TRUE, FALSE))
operation <- operation * Ovariable(output = data.frame(Plant = "Operation", Result = 1), marginal = c(TRUE, FALSE))
cost <- combine(investment, maintenance, operation)
cost <- OpasnetUtils::combine(investment, maintenance, operation)
marginals <- character()
marginals <- character()
for(i in colnames(cost@output)[cost@marginal]) {
for(i in colnames(cost@output)[cost@marginal]) {
Line 869: Line 536:
</rcode>
</rcode>


==== Original version ====
=== Fuel use and fuel shares in generic processes ===
 
There is an alternative way for calculating fuel use. It is based on the idea that for each heating type, there is a constant share of fuels used. For some heating types, this is generic and is shown on this page. For some others, the constant is case-specific and is determined on a case-specific page.
 
The table below contains connections of heating types and fuel usage in generic situations. There may be case-specific differences, which must be handled separately.
 
<t2b name='Fuel use in different heating types' index='Heating,Burner,Fuel' obs='Fraction' desc='Description' unit='-'>
Wood|Household|Wood|1|
Oil|Household|Light oil|1|
Gas|Household|Gas|1|
Heating oil|Household|Light oil|1|
Coal|Household|Coal|1|
Other sources|Household|Other sources|1|
No energy source|Household|Other sources|1|
Geothermal|Grid|Electricity|0.3|Geothermal does not sum up to 1 because more heat is produced than electricity consumed.
Centrifuge, hydro-extractor|Grid|Electricity|0.3|Not quite clear what this is but presumably a heat pump.
Solar heater/ collector|Grid|Electricity|0.1|Use only; life-cycle impacts omitted.
Electricity|Grid|Electricity|1|
District|Undefined|Heat|1|
</t2b>
 
<rcode name='fuelSharesgeneric' embed=1 label='Initiate fuelSharesgeneric (only for developers)'>
 
# This is code Op_en5141/fuelSharesgeneric (only generic) on page [[Energy balance]].
library(OpasnetUtils)
 
fuelSharesgeneric <- Ovariable("fuelSharesgeneric", ddata = "Op_en5141", subset = "Fuel use in different heating types") # [[Energy balance]]
colnames(fuelSharesgeneric@data) <- gsub("[ \\.]", "_", colnames(fuelSharesgeneric@data))
#fuelSharesgeneric@data <- merge(fuelSharesgeneric@data, data.frame(Time = 1900:2080))
 
objects.store(fuelSharesgeneric)
cat("Object fuelSharesgeneric initiated!\n")
 
</rcode>
 
<rcode name="fuelUse" label="Initiate fuelUse (developers only)" embed=1>
## This is code Op_en5141/fuelUse on page [[Energy balance]].
library(OpasnetUtils)
 
fuelUse <- Ovariable("fuelUse",
dependencies = data.frame(
Name = c("energyUse", "fuelShares")
),
formula = function(...) {
out <- energyUse * fuelShares
return(out)
}
)
 
objects.store(fuelUse)
cat("Ovariable fuelUse stored.\n")
</rcode>
<noinclude>
 
=== Old version with a set of linear equations ===
 
* Energy balances are described as input = output on a coarse level (called classes) where the structure is the same or similar to the OECD energy balance tables. If possible, this is described on the Energy balance method level and it is shared by all cities.
* On more detailed (variable level in the matrix), the fraction of each variable of the total class are described separately. Fractions are city specific and they are described on city level in a separate table.
* Based on the fraction table, detailed equations with variables are created. The format will be fraction * class total = variable.
* The last fraction has zero degrees of freedom when the class total is given. However, it must have a variable and thus a row in the fraction table. The result for that variable is an empty cell (which results in NA).
* Unlike in the previous version, all variables are given either as values or equations, and the user interface is not used for BAU. In contrast, user interface or decision table may be used to derive values for alternative scenarios.
* To make this work, the city-specific fraction data must be defined as ovariable (so that it can be changed with a decision table), and also the energy balance method must be described asa ovariable. How are we going to make the two interplay, as we may want to have several cities?
** Define one city ovariable and evaluate energy balance with that. The ovariable has a generic name. Then, define a new city ovariable with the same name and re-evaluate the energy balance ovariable; this must be done so that the two cities are appended rather than replaced.
** city ovariables are appended first into a large fraction table, and then that is used to create the large energy balance matrix. {{defend|# |This is clearly better.|--[[User:Jouni|Jouni]] 17:09, 21 February 2013 (EET)}}
* The city-specific ovariable may have Iter and other indices. A separate matrix is created and solved for each unique combination of indices. This makes it possible to have a very flexible approach.
* We should check if the energy balance matrix (see Matti's Excel) has city-specific equations. If possible, energy transformations are described as generic equations on the energy balance method.
* Structure of OECD Energy balance tables (data):
** Fuel (given as observation columns in OECD table)
** Activity (row in OECD table)
** Description
* Structure of the generic process table
** Equation,
** Col,
** Result,
** Description? {{attack|# |This does not join up in a coherent way.|--[[User:Jouni|Jouni]] 17:09, 21 February 2013 (EET)}}
* Columns for fraction table
** Class
** Item
** Result (fraction)
** Indices as needed
 
; Example table for making matrices from text format equations. CHPcapacity describes which of the piecewise linear equations should be used. Policy is a decision option that alters the outcome. Dummy is only for compatibility but it is not used.
 
<t2b name="Equations" index="CHPcapacity,Policy,Equation" obs="Dummy" desc="Description" unit="GWh /a">
|Biofuel|CHP renewable = CHP peat|1|Biofuel policy contains half biofuels, half peat
|BAU|CHP renewable = 89.24|1|
||CHP peat + CHP renewable + CHP oil = CHP heat + CHP electricity|1|
||CHP peat = 90-98*CHP oil|1|
||CHP electricity = 0.689*CHP heat|1|
CHP<1000||H heat = 0.08*CHP heat|1|Small heat plants reflect the total heat need
CHP>1000||CHP heat + CHP electricity = 1000|1|But production capacity of CHP may be overwhelmed, decoupling CHP heat and H heat.
||H biogas + H oil = H heat|1|
||H oil = 18.973*H biogas|1|
||Bought electricity + CHP electricity = Cons electricity|1|
||CHP heat + H heat = Cons heat|1|
||Cons electricity = 900-1100|1|
||Cons heat = 900-1000|1|
</t2b>
 
; Example table to describe the details about nonlinear equations.
 
<t2b name="Nonlinearity parameters" index="critVar,critIndex,rescol,critLocLow,critLocHigh" obs="critValue" unit="GWh /a">
Cons heat|CHPcapacity|Result|CHP<1000|CHP>1000|1080
</t2b>
 
 
; This table is fetched if there are no nonlinearities. Therefore, there is no need to copy it to the case study page.
 
<t2b name="No nonlinearities" index="critVar,critIndex,rescol,critLocLow,critLocHigh" obs="critValue" unit="GWh /a">
|||||
</t2b>
 
; This table is fetched if there are no modelled upstream variables that would affect the equations.
 
<t2b name="No modelled upstream variables" index="energybalanceVars" obs="Result" unit="-">
|
</t2b>
 


Stored objects below used by [[Energy balance in Kuopio]].
Stored objects below used by [[Energy balance in Kuopio]].
Line 1,089: Line 873:
{{Helsinki energy decision 2015}}
{{Helsinki energy decision 2015}}


* A previous model version written by Jouni. It used linear optimising but was not fully developed and was not actually used in any final assessment. [http://en.opasnet.org/en-opwiki/index.php?title=Energy_balance&oldid=38536#Calculations]
* A new version of energy balance is using [https://stat.ethz.ch/R-manual/R-devel/library/stats/html/optim.html optim function]. See [[Energy balance in Helsinki]].
* A new version of energy balance is using [https://stat.ethz.ch/R-manual/R-devel/library/stats/html/optim.html optim function]. See [[Energy balance in Helsinki]].
* [http://www.nzz.ch/meinung/kommentare/energiepolitik-bitte-wenden-1.18451262 Energiepolitik, bitte wenden!] Neue Zürcher Zeitung 27.12.2014 by Giorgio V. Müller
* [http://www.nzz.ch/meinung/kommentare/energiepolitik-bitte-wenden-1.18451262 Energiepolitik, bitte wenden!] Neue Zürcher Zeitung 27.12.2014 by Giorgio V. Müller
Line 1,120: Line 905:


==Related files==
==Related files==
</noinclude>

Latest revision as of 07:49, 10 April 2019



Question

What is energy balance and how is it modelled?

Answer

Summing up the amount of energy produced and subtracting the amount of energy consumed within a time period gives the energy balance. Since the electricity grid and district heat network lack significant storage mechanics, the balance has to be virtually zero over short periods. When considering the balance of a particular area (e.g. Helsinki), we can make the assumption that electricity can be imported and exported in international markets. The energy in the district heat network, however, has to be produced locally. This sets up the non-trivial problem of optimising production so that there are no significant deficits as well as minimising losses and maximising profits. This problem is solved (to some extent) by market forces in the real world.

In Opasnet, there are two different ways to calculate energy balance. Our most recent energy balance model uses linear programming tools to solve an optimum for the activity of a given set of production units in simulated instances created by the main model. The main model is responsible for the decision making aspects, while the energy balance optimisation only functions as an approximation of real world market mechanics. This version was used in Helsinki energy decision 2015.

The previous version was based on setting up a set of linear equations describing the inputs, outputs, and shares of different energy and plant processes. This approach is less flexible, because it does not use an optimising function and everything must be described as linear (or piecewise linear). However, this approach was successfully used in Energy balance in Kuopio and Energy balance in Suzhou.

This code is an example how the energy balance model is used in a city case. The data comes from Helsinki energy decision 2015.

+ Show code

Rationale

Energy balance with linear programming

The linear programming problem is set up as follows.

For each production unit: let xi be activity of the plant. Lets also have variables yj for deficits and excesses for each type of energy produced.

The objective function is the function we are optimising. Each production unit has a unit profit per activity denoted by ai which is determined by the amount of different input commodities (e.g. coal) per amount of different output commodities (i.e. electricity and heat) and their market prices. Also, lets say we want to make sure that district heat demand is always met when possible and have a large penalty factor for each unit of heat demand not met (1 M€ in the model). In addition, it must be noted that excess district heat becomes wasted so it counts as loss. Let these deficit and excess related losses be denoted by bj. The whole objective function then becomes: sum(xiai) + sum(yjbj).

The values of variables are constrained by equalities and inequalities: the sum of production of a commodity is equal to its demand minus deficit plus excess, activity is constrained by the maximum capacity and all variables are non-negative by definition. This can be efficiently solved by computers for each given instance. Production wind-up and wind-down is ignored, since time continuity is not considered. As a consequence fuel limits (e.g. diminishing hydropower capacity) are not modelled completely either.

Ovariables like EnergyNetworkOptim below are used in Helsinki energy decision 2015. Prices of fuels in heat production are used as direct inputs in the optimising.

+ Show code

Fuel use and fuel shares in generic processes

There is an alternative way for calculating fuel use. It is based on the idea that for each heating type, there is a constant share of fuels used. For some heating types, this is generic and is shown on this page. For some others, the constant is case-specific and is determined on a case-specific page.

The table below contains connections of heating types and fuel usage in generic situations. There may be case-specific differences, which must be handled separately.

Fuel use in different heating types(-)
ObsHeatingBurnerFuelFractionDescription
1WoodHouseholdWood1
2OilHouseholdLight oil1
3GasHouseholdGas1
4Heating oilHouseholdLight oil1
5CoalHouseholdCoal1
6Other sourcesHouseholdOther sources1
7No energy sourceHouseholdOther sources1
8GeothermalGridElectricity0.3Geothermal does not sum up to 1 because more heat is produced than electricity consumed.
9Centrifuge, hydro-extractorGridElectricity0.3Not quite clear what this is but presumably a heat pump.
10Solar heater/ collectorGridElectricity0.1Use only; life-cycle impacts omitted.
11ElectricityGridElectricity1
12DistrictUndefinedHeat1

+ Show code

+ Show code


Old version with a set of linear equations

  • Energy balances are described as input = output on a coarse level (called classes) where the structure is the same or similar to the OECD energy balance tables. If possible, this is described on the Energy balance method level and it is shared by all cities.
  • On more detailed (variable level in the matrix), the fraction of each variable of the total class are described separately. Fractions are city specific and they are described on city level in a separate table.
  • Based on the fraction table, detailed equations with variables are created. The format will be fraction * class total = variable.
  • The last fraction has zero degrees of freedom when the class total is given. However, it must have a variable and thus a row in the fraction table. The result for that variable is an empty cell (which results in NA).
  • Unlike in the previous version, all variables are given either as values or equations, and the user interface is not used for BAU. In contrast, user interface or decision table may be used to derive values for alternative scenarios.
  • To make this work, the city-specific fraction data must be defined as ovariable (so that it can be changed with a decision table), and also the energy balance method must be described asa ovariable. How are we going to make the two interplay, as we may want to have several cities?
    • Define one city ovariable and evaluate energy balance with that. The ovariable has a generic name. Then, define a new city ovariable with the same name and re-evaluate the energy balance ovariable; this must be done so that the two cities are appended rather than replaced.
    • city ovariables are appended first into a large fraction table, and then that is used to create the large energy balance matrix. ←--#: . This is clearly better. --Jouni 17:09, 21 February 2013 (EET) (type: truth; paradigms: science: defence)
  • The city-specific ovariable may have Iter and other indices. A separate matrix is created and solved for each unique combination of indices. This makes it possible to have a very flexible approach.
  • We should check if the energy balance matrix (see Matti's Excel) has city-specific equations. If possible, energy transformations are described as generic equations on the energy balance method.
  • Structure of OECD Energy balance tables (data):
    • Fuel (given as observation columns in OECD table)
    • Activity (row in OECD table)
    • Description
  • Structure of the generic process table
    • Equation,
    • Col,
    • Result,
    • Description? ⇤--#: . This does not join up in a coherent way. --Jouni 17:09, 21 February 2013 (EET) (type: truth; paradigms: science: attack)
  • Columns for fraction table
    • Class
    • Item
    • Result (fraction)
    • Indices as needed
Example table for making matrices from text format equations. CHPcapacity describes which of the piecewise linear equations should be used. Policy is a decision option that alters the outcome. Dummy is only for compatibility but it is not used.
Equations(GWh /a)
ObsCHPcapacityPolicyEquationDummyDescription
1BiofuelCHP renewable = CHP peat1Biofuel policy contains half biofuels, half peat
2BAUCHP renewable = 89.241
3CHP peat + CHP renewable + CHP oil = CHP heat + CHP electricity1
4CHP peat = 90-98*CHP oil1
5CHP electricity = 0.689*CHP heat1
6CHP<1000H heat = 0.08*CHP heat1Small heat plants reflect the total heat need
7CHP>1000CHP heat + CHP electricity = 10001But production capacity of CHP may be overwhelmed, decoupling CHP heat and H heat.
8H biogas + H oil = H heat1
9H oil = 18.973*H biogas1
10Bought electricity + CHP electricity = Cons electricity1
11CHP heat + H heat = Cons heat1
12Cons electricity = 900-11001
13Cons heat = 900-10001
Example table to describe the details about nonlinear equations.
Nonlinearity parameters(GWh /a)
ObscritVarcritIndexrescolcritLocLowcritLocHighcritValue
1Cons heatCHPcapacityResultCHP<1000CHP>10001080


This table is fetched if there are no nonlinearities. Therefore, there is no need to copy it to the case study page.
No nonlinearities(GWh /a)
ObscritVarcritIndexrescolcritLocLowcritLocHighcritValue
1
This table is fetched if there are no modelled upstream variables that would affect the equations.
No modelled upstream variables(-)
ObsenergybalanceVarsResult
1


Stored objects below used by Energy balance in Kuopio.

+ Show code

How to give uncertain parameters?

  • In equations, the content is interpreted only inside solveMatrix. Therefore, the typical approach where all unique index combinations are run one at a time does not work.
  • There should be an update in parameter interpretation for terms with one entry only. It can no longer be based on as.numeric, if distributions (=text) is allowed.
    • If it starts with [a-z.] it is a variable name.
    • If it starts with [0-9<\\-] it is a parameter value.
  • Instead of params[[i]] and [[vars]] vectors, a data.frame will be created with Result as the params column.
  • The data.frame is then interpreted with N = N. If parameters are probabilistic, Iter column will appear.
  • When all parameters have been interpreted, check if Iter exists.
  • If Iter exists, make a for loop for all values of Iter.
    • Create a matrix from the parameters and solve.
    • Rbind the result to a data.frame with Iter.
  • Return the output.
  • Old code with an input table with columns Equation, Col, Result, Description: [1]

See also

Helsinki energy decision 2015
In English
Assessment Main page | Helsinki energy decision options 2015
Helsinki data Building stock in Helsinki | Helsinki energy production | Helsinki energy consumption | Energy use of buildings | Emission factors for burning processes | Prices of fuels in heat production | External cost
Models Building model | Energy balance | Health impact assessment | Economic impacts
Related assessments Climate change policies in Helsinki | Climate change policies and health in Kuopio | Climate change policies in Basel
In Finnish
Yhteenveto Helsingin energiapäätös 2015 | Helsingin energiapäätöksen vaihtoehdot 2015 | Helsingin energiapäätökseen liittyviä arvoja | Helsingin energiapäätös 2015.pptx
Urgenche research project 2011 - 2014: city-level climate change mitigation
Urgenche pages

Urgenche main page · Category:Urgenche · Urgenche project page (password-protected)

Relevant data
Building stock data in Urgenche‎ · Building regulations in Finland · Concentration-response to PM2.5 · Emission factors for burning processes · ERF of indoor dampness on respiratory health effects · ERF of several environmental pollutions · General criteria for land use · Indoor environment quality (IEQ) factors · Intake fractions of PM · Land use in Urgenche · Land use and boundary in Urgenche · Energy use of buildings

Relevant methods
Building model · Energy balance · Health impact assessment · Opasnet map · Help:Drawing graphs · OpasnetUtils‎ · Recommended R functions‎ · Using summary tables‎

City Kuopio
Climate change policies and health in Kuopio (assessment) · Climate change policies in Kuopio (plausible city-level climate policies) · Health impacts of energy consumption in Kuopio · Building stock in Kuopio · Cost curves for energy (prioritization of options) · Energy balance in Kuopio (energy data) · Energy consumption and GHG emissions in Kuopio by sector · Energy consumption classes (categorisation) · Energy consumption of heating of buildings in Kuopio · Energy transformations (energy production and use processes) · Fuels used by Haapaniemi energy plant · Greenhouse gas emissions in Kuopio · Haapaniemi energy plant in Kuopio · Land use in Kuopio · Building data availability in Kuopio · Password-protected pages: File:Heat use in Kuopio.csv · Kuopio housing

City Basel
Buildings in Basel (password-protected)

Energy balances
Energy balance in Basel · Energy balance in Kuopio · Energy balance in Stuttgart · Energy balance in Suzhou


References


Related files