From a6d10650705d7b0e110ffb5f1a3b09f735f5e371 Mon Sep 17 00:00:00 2001
From: R0slyn <roslyn.henry.08@aberdeen.ac.uk>
Date: Tue, 17 Oct 2017 16:22:52 +0100
Subject: [PATCH] dietary_branch
---
GAMS/IntExtOpt.gms | 7 +-
src/ac/ed/lurg/ModelConfig.java | 12 +++-
src/ac/ed/lurg/ModelMain.java | 46 +++++++++++-
src/ac/ed/lurg/country/CountryAgent.java | 24 +++++++
.../country/gams/GamsLocationOptimiser.java | 6 +-
.../ed/lurg/demand/AbstractDemandManager.java | 72 ++++++++++++++++---
src/ac/ed/lurg/types/CommodityType.java | 26 +++++--
7 files changed, 170 insertions(+), 23 deletions(-)
diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms
index e1b531c6..9a011219 100644
--- a/GAMS/IntExtOpt.gms
+++ b/GAMS/IntExtOpt.gms
@@ -231,7 +231,7 @@ $gdxin
parameter totalCropland(location);
parameter netImportAmount(all_types);
parameter netImportCost(all_types);
- parameter feedCost(all_types);
+ parameter feedCostRate(feed_crop);
* Production quantities based on smaller area (before unhandledCropArea adjustment applied)
totalProd(crop) = sum(location, area.l(crop, location) * yield.l(crop, location));
@@ -244,9 +244,14 @@ $gdxin
totalCropland(location) = sum(crop_less_pasture, area.l(crop_less_pasture, location));
totalArea(crop) = sum(location, area.l(crop, location));
+ feedCostRate(feed_crop)$[totalProd(feed_crop) <> 0] = totalProdCost(feed_crop) / totalProd(feed_crop);
+ totalProdCost('ruminants') = sum(feed_crop, ruminantFeed.l(feed_crop) * feedCostRate(feed_crop));
+ totalProdCost('monogastrics') = sum(feed_crop, monogastricFeed.l(feed_crop) * feedCostRate(feed_crop));
+
netImportCost(import_crop) = importAmount.l(import_crop) * importPrices(import_crop) - exportAmount.l(import_crop) * exportPrices(import_crop);
netImportAmount(import_crop) = importAmount.l(import_crop) - exportAmount.l(import_crop);
Scalar ms 'model status', ss 'solve status';
ms=LAND_USE.modelstat;
ss=LAND_USE.solvestat;
+
\ No newline at end of file
diff --git a/src/ac/ed/lurg/ModelConfig.java b/src/ac/ed/lurg/ModelConfig.java
index 14713281..5d4b01de 100644
--- a/src/ac/ed/lurg/ModelConfig.java
+++ b/src/ac/ed/lurg/ModelConfig.java
@@ -163,6 +163,7 @@ public class ModelConfig {
public static final String LAND_COVER_OUTPUT_FILE = OUTPUT_DIR + File.separator + "lc.txt";
public static final String PRICES_OUTPUT_FILE = OUTPUT_DIR + File.separator + "prices.txt";
public static final String DEMAND_OUTPUT_FILE = OUTPUT_DIR + File.separator + "demand.txt";
+ public static final String DOMESTIC_OUTPUT_FILE = OUTPUT_DIR + File.separator + "domestic.txt";
public static final boolean OUTPUT_FOR_LPJG = getBooleanProperty("OUTPUT_FOR_LPJG", true);
public static final boolean INTERPOLATE_OUTPUT_YEARS = getBooleanProperty("INTERPOLATE_OUTPUT_YEARS", true);
@@ -238,7 +239,14 @@ public class ModelConfig {
public static final String BIOENERGY_DEMAND_SCENARIO = getProperty("BIOENERGY_DEMAND_SCENARIO", "BAU");
public static final double BIOENERGY_DEMAND_SHIFT = IS_CALIBRATION_RUN ? 1.0 : getDoubleProperty("BIOENERGY_DEMAND_SHIFT", 1.0);
// public static final double BIOENERGY_HEATING_VALUE_GJ_PER_T = getDoubleProperty("BIOENERGY_HEATING_VALUE_GJ_PER_T", 17.5); // GJ per t DM
-
+
+ //Dietary stuff
+ public static final boolean ACTIVE_DIETARY_SHIFT = getBooleanProperty("ACTIVE_DIETARY_SHIFT", false);
+ public static final int DIETARY_CHANGE_START_YEAR = getIntProperty("DIETARY_CHANGE_START_YEAR", 2010);
+ public static final int DIETARY_CHANGE_END_YEAR = getIntProperty("DIETARY_CHANGE_END_YEAR", 2050);
+ public static final double DIETARY_CHANGE_ANNUAL_RATE = getDoubleProperty("DIETARY_CHANGE_ANNUAL_RATE", 0.05);
+
+
public static final double MARKET_LAMBA = getDoubleProperty("MARKET_LAMBA", 0.4); // controls international market price adjustment rate
public static final double POPULATION_AGGREG_LIMIT = getDoubleProperty("POPULATION_AGGREG_LIMIT", 30.0); // in millions, smaller countries are aggregated on a regional basis
@@ -274,4 +282,6 @@ public class ModelConfig {
public static final int NUM_YIELD_CLUSTERS = getIntProperty("NUM_YIELD_CLUSTERS", 5000);
public static final long RANDOM_SEED = getIntProperty("RANDOM_SEED", 1974329); // any number will do
+
+
}
diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java
index 37dbdd75..05206b69 100644
--- a/src/ac/ed/lurg/ModelMain.java
+++ b/src/ac/ed/lurg/ModelMain.java
@@ -150,7 +150,7 @@ public class ModelMain {
double gen2Increase = (gen2EcDDemand>previousGen2EcDDemand) ? gen2EcDDemand - previousGen2EcDDemand : 0.0;
CropToDoubleMap totalImportCommodities = new CropToDoubleMap();
CropToDoubleMap totalExportCommodities = new CropToDoubleMap();
-
+
for (CountryAgent ca : countryAgents) {
LogWriter.println("Country " + ca.getCountry());
@@ -206,6 +206,7 @@ public class ModelMain {
else
totalExportCommodities.incrementValue(c, -countryNetImports);
}
+
}
// energycrops are a special case where demand in global and exogenously specified
@@ -324,12 +325,51 @@ public class ModelMain {
LogWriter.print(e);
}
}
+
+
+ private void writeDomesticProductionFile(Timestep timestep) {
+ try {
+ StringBuffer sbHeadings = new StringBuffer("Year, Country, Crop,Production cost, export price, import price");
+ BufferedWriter outputFile = getFileWriter(timestep, ModelConfig.DOMESTIC_OUTPUT_FILE,
+ sbHeadings.toString());
+
+ for (CropType crop : CropType.getAllItems()) {
+ for (CountryAgent country : countryAgents) {
+ LogWriter.println("Country " + country.getCountry().getName() + " crop " + crop);
+ Double importPrice;
+ Double exportPrice;
+ Double prodCosts = country.getDomesticProductionCosts().get(crop);
+
+ if (crop.isImportedCrop()) {
+ importPrice = country.getCurrentCountryPrices().get(crop).getImportPrice();
+ exportPrice = country.getCurrentCountryPrices().get(crop).getExportPrice();
+ } else {
+ importPrice = null;
+ exportPrice = null;
+ }
+ StringBuffer sbData = new StringBuffer();
+ sbData.append(String.format("%d,%s, %s", timestep.getYear(), country.getCountry(),crop.getGamsName()));
+ sbData.append(String.format(",%.3f,%.3f,%.3f", prodCosts, importPrice, exportPrice));
+
+ outputFile.write(sbData.toString());
+ outputFile.newLine();
+ }
+
+ }
+ outputFile.close();
+
+ } catch (IOException e) {
+ LogWriter.print(e);
+ }
+
+ }
private void outputTimestepResults(Timestep timestep, RasterSet<LandUseItem> landUseRaster, YieldRaster yieldSurfaces) {
writeLandCoverFile(timestep, landUseRaster);
writeGlobalMarketFile(timestep);
writeDemandFile(timestep);
+ writeDomesticProductionFile(timestep);
if (ModelConfig.OUTPUT_FOR_LPJG) {
for (int outputYear : timestep.getYearsFromLast()) {
@@ -374,6 +414,8 @@ public class ModelMain {
// outputLandCover(timestep.getYear(), landUseRaster, LandCoverType.PASTURE);
}
+
+
private void outputWaterAvailablity(Timestep timestep, IrrigationRasterSet irrigiationRS) {
new RasterOutputer<Double, IrrigationItem>(irrigiationRS, ModelConfig.OUTPUT_DIR + File.separator + timestep.getYear() + File.separator + "IrrigConstraint.asc") {
@Override
@@ -528,7 +570,6 @@ public class ModelMain {
areasItem.setProtectedArea(landCover.getProtectedArea());
landUseRaster.put(key, areasItem);
}
-
return landUseRaster;
}
@@ -566,4 +607,5 @@ public class ModelMain {
initialStockLevels.put(CropType.ENERGY_CROPS, 0.0);
return initialStockLevels;
}
+
}
diff --git a/src/ac/ed/lurg/country/CountryAgent.java b/src/ac/ed/lurg/country/CountryAgent.java
index abaa9f9a..12292c04 100644
--- a/src/ac/ed/lurg/country/CountryAgent.java
+++ b/src/ac/ed/lurg/country/CountryAgent.java
@@ -41,6 +41,7 @@ public class CountryAgent {
private Map<CropType, CountryPrice> currentCountryPrices;
private Map<CropType, Double> tradeBarriers;
private RasterSet<IntegerRasterItem> yieldClusters;
+ private Map<CropType, Double> domesticProductionCosts;
public CountryAgent(AbstractDemandManager demandManager,CompositeCountry country, RasterSet<LandUseItem> cropAreaRaster,
Map<CropType, CropUsageData> cropUsageData, Map<CropType, Double> tradeBarriers, RasterSet<IntegerRasterItem> yieldClusters) {
@@ -52,6 +53,7 @@ public class CountryAgent {
GamsRasterOutput initialData = new GamsRasterOutput(cropAreaRaster, cropUsageData);
previousGamsRasterOutput = initialData;
+
}
public CompositeCountry getCountry() {
@@ -125,7 +127,10 @@ public class CountryAgent {
LogWriter.println("Running " + country.getName() + ", currentTimestep " + currentTimestep);
GamsRasterOutput result = opti.run();
+ domesticProductionCosts = calculateDomesticProductionCosts(result);
+
previousGamsRasterOutput = result;
+
return result;
}
@@ -135,6 +140,10 @@ public class CountryAgent {
public Map<CommodityType, Double> getCurrentProjectedDemand() {
return currentProjectedDemand;
}
+
+ public Map<CropType, CountryPrice> getCurrentCountryPrices() {
+ return currentCountryPrices;
+ }
private GamsRasterInput getGamsRasterInput(RasterSet<IrrigationItem> irrigData, YieldRaster countryYieldSurfaces, double gen2EcIncrease) {
double allowedImportChange;
@@ -191,4 +200,19 @@ public class CountryAgent {
return countryPrices;
}
+
+ public Map<CropType, Double> getDomesticProductionCosts() {
+ return domesticProductionCosts;
+ }
+
+ private Map<CropType, Double> calculateDomesticProductionCosts(GamsRasterOutput result) {
+
+ Map<CropType, Double> domesticProdCosts = new HashMap<CropType, Double>();
+ Map <CropType, CropUsageData> cropUsage = result.getCropUsageData();
+
+ for (Entry<CropType, CropUsageData> entry : cropUsage.entrySet()) {
+ domesticProdCosts.put(entry.getKey(),entry.getValue().getProdCost());
+ }
+ return domesticProdCosts;
+ }
}
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
index fc385087..515ee338 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
@@ -364,10 +364,10 @@ public class GamsLocationOptimiser {
for (CropType meatTypes : CropType.getMeatTypes()) {
netImport = getParmValue(parmNetImports, meatTypes.getGamsName());
prod = getParmValue(parmProd, meatTypes.getGamsName());
-// prodCost = getParmValue(parmProdCost, meatTypes.getGamsName()); currently not calculated in GAMS correctly, not sure it's needed
+ prodCost = getParmValue(parmProdCost, meatTypes.getGamsName()); // currently not calculated in GAMS correctly, not sure it's needed
- cropUsageData.put(meatTypes, new CropUsageData(0.0, 0.0, netImport, prod, Double.NaN));
- if (DEBUG) LogWriter.println(String.format("\n%s:\t\t\t\t\tnetImports= %.3f,\tprod= %.3f", meatTypes.getGamsName(), netImport, prod));
+ cropUsageData.put(meatTypes, new CropUsageData(0.0, 0.0, netImport, prod, prodCost));
+ if (DEBUG) LogWriter.println(String.format("\n%s:\t\t\t\t\tnetImports= %.3f,\tprod= %.3f,\tprodCost= %.3f", meatTypes.getGamsName(), netImport, prod, prodCost));
}
LogWriter.println(String.format("\n%s %s: Total area= %.1f (crop=%.1f, pasture %.1f)",
inputData.getCountryInput().getCountry(), inputData.getTimestep().getYear(),
diff --git a/src/ac/ed/lurg/demand/AbstractDemandManager.java b/src/ac/ed/lurg/demand/AbstractDemandManager.java
index d18383ee..4efe551d 100644
--- a/src/ac/ed/lurg/demand/AbstractDemandManager.java
+++ b/src/ac/ed/lurg/demand/AbstractDemandManager.java
@@ -7,8 +7,10 @@ import ac.ed.lurg.ModelConfig;
import ac.ed.lurg.Timestep;
import ac.ed.lurg.country.CompositeCountry;
import ac.ed.lurg.country.CompositeCountryManager;
+import ac.ed.lurg.country.GlobalPrice;
import ac.ed.lurg.country.SingleCountry;
import ac.ed.lurg.types.CommodityType;
+import ac.ed.lurg.types.CropType;
import ac.ed.lurg.utils.LogWriter;
public abstract class AbstractDemandManager {
@@ -18,37 +20,45 @@ public abstract class AbstractDemandManager {
public AbstractDemandManager(CompositeCountryManager compositeCountryManager) {
this.compositeCountryManager = compositeCountryManager;
- bioenergyDemandManager = new BioenergyDemandManager();
+ bioenergyDemandManager = new BioenergyDemandManager();
}
public Map<CommodityType, Double> getDemand(CompositeCountry cc, int year) {
- if (!ModelConfig.CHANGE_DEMAND_YEAR)
+ if (!ModelConfig.CHANGE_DEMAND_YEAR)
year = ModelConfig.BASE_YEAR;
Map<CommodityType, Double> demandMap = new HashMap<CommodityType, Double>();
-
+ Map<CommodityType, Double> foodDemandMap;
+
for (SingleCountry c : compositeCountryManager.getAllForCompositeCountry(cc)) {
- Map<CommodityType, Double> foodDemandMap = getFoodDemand(c, year);
+
+ if (ModelConfig.ACTIVE_DIETARY_SHIFT) {
+ foodDemandMap = calculateDietaryShift(c, year);
+ } else
+ foodDemandMap = getFoodDemand(c, year);
for (CommodityType commodity : CommodityType.getAllItems()) {
+
Double food = foodDemandMap.get(commodity);
double bioenergy = bioenergyDemandManager.getFirstGenBioenergyDemand(c, year, commodity);
- // LogWriter.println(String.format("Country %s comm %s: %f", country, commodity, bioenergy));
+ // LogWriter.println(String.format("Country %s comm %s: %f",
+ // country, commodity, bioenergy));
- double newDemand = ((food==null ? 0.0 : food.doubleValue()) + bioenergy) * (1.0 + commodity.getSeedWasteRate());
+ double newDemand = ((food == null ? 0.0 : food.doubleValue()) + bioenergy)
+ * (1.0 + commodity.getSeedWasteRate());
Double prevDemand = demandMap.get(commodity);
if (prevDemand != null)
- newDemand += prevDemand.doubleValue(); // aggregate if composite country
+ newDemand += prevDemand.doubleValue(); // aggregate if
+ // composite country
demandMap.put(commodity, newDemand);
}
}
-
return demandMap;
}
-
+
abstract Map<CommodityType, Double> getFoodDemand(SingleCountry c, int year);
public double getSecondGenBioenergyDemand(Timestep timestep) {
@@ -57,4 +67,48 @@ public abstract class AbstractDemandManager {
LogWriter.println("Global gen2 bioenergy demand in " + timestep + " is " + d);
return d;
}
+
+ private Map<CommodityType, Double> calculateDietaryShift(SingleCountry c, int year) {
+
+ int yearsOfChange = Math.min(ModelConfig.DIETARY_CHANGE_END_YEAR - ModelConfig.DIETARY_CHANGE_START_YEAR,
+ year - ModelConfig.DIETARY_CHANGE_START_YEAR);
+ double dietaryDemandAdj = yearsOfChange > 0 ? (yearsOfChange * ModelConfig.DIETARY_CHANGE_ANNUAL_RATE) : 0.0;
+
+ Map<CommodityType, Double> originalFoodDemandMap = getFoodDemand(c, year);
+
+ Double monogastricDemand = originalFoodDemandMap.get(CommodityType.MONOGASTRICS);
+ double mongoastricCalorieShortfall = (monogastricDemand != null)
+ ? monogastricDemand * dietaryDemandAdj * CommodityType.MONOGASTRICS.getEnergyPerTon() : 0.0;
+
+ Double ruminantDemand = originalFoodDemandMap.get(CommodityType.RUMINANTS);
+ double ruminantsCalorieShortfall = (ruminantDemand != null)
+ ? ruminantDemand * dietaryDemandAdj * CommodityType.RUMINANTS.getEnergyPerTon() : 0.0;
+
+ Map<CommodityType, Double> updatedFoodDemandMap = new HashMap<CommodityType, Double>();
+
+ for (CommodityType commodity : CommodityType.getAllItems()) {
+ Double originalCommodityDemand = originalFoodDemandMap.get(commodity);
+
+ if (originalCommodityDemand != null) {
+ if (commodity.isAnimalProduct()) {
+ updatedFoodDemandMap.put(commodity,
+ originalCommodityDemand - originalCommodityDemand * dietaryDemandAdj);
+ } else {
+ updatedFoodDemandMap.put(commodity,
+ originalCommodityDemand
+ + commodity.getMeatSubstitutionProportion() * ruminantsCalorieShortfall
+ / commodity.getEnergyPerTon()
+ + commodity.getMeatSubstitutionProportion() * mongoastricCalorieShortfall
+ / commodity.getEnergyPerTon());
+
+ }
+ }
+
+ }
+
+// LogWriter.println("original map " + originalFoodDemandMap);
+// LogWriter.println("updated map " + updatedFoodDemandMap);
+ return updatedFoodDemandMap;
+ }
+
}
diff --git a/src/ac/ed/lurg/types/CommodityType.java b/src/ac/ed/lurg/types/CommodityType.java
index 615b6efb..b66ff622 100644
--- a/src/ac/ed/lurg/types/CommodityType.java
+++ b/src/ac/ed/lurg/types/CommodityType.java
@@ -8,21 +8,25 @@ import ac.ed.lurg.ModelConfig;
public enum CommodityType {
- CEREALS("Cereals", "cereals", false),
- OILCROPS("Oilcrops", "oilcrops", false),
- PULSES("Pulses", "pulses", false),
- STARCHY_ROOTS("Starchy Roots", "starchyRoots", false),
- MONOGASTRICS("Monogastrics", "monogastrics", true),
- RUMINANTS("Ruminants", "ruminants", true);
+ CEREALS("Cereals", "cereals", false, 0.9393854,0.2),
+ OILCROPS("Oilcrops", "oilcrops", false, 2.0640544,0.1),
+ PULSES("Pulses", "pulses", false, 2.5886803,0.4),
+ STARCHY_ROOTS("Starchy Roots", "starchyRoots", false, 0.4925546,0.3),
+ MONOGASTRICS("Monogastrics", "monogastrics", true, 0.4680600, 1.0),
+ RUMINANTS("Ruminants", "ruminants", true, 0.6135583, 1.0);
private String faoName;
private String gamsName;
private boolean isAnimalProduct;
+ private double energyPerTon;
+ private double meatSubstitutionProportion;
- CommodityType (String faoName, String gamsName, boolean isAnimalProduct) {
+ CommodityType (String faoName, String gamsName, boolean isAnimalProduct, double energyPerTon, double meatSubstitutionProportion) {
this.faoName = faoName;
this.gamsName = gamsName;
this.isAnimalProduct = isAnimalProduct;
+ this.energyPerTon = energyPerTon;
+ this.meatSubstitutionProportion = meatSubstitutionProportion;
}
private static final Map<String, CommodityType> faoNameCache = new HashMap<String, CommodityType>();
@@ -57,6 +61,14 @@ public enum CommodityType {
return isAnimalProduct;
}
+ public double getEnergyPerTon(){
+ return energyPerTon;
+ }
+
+ public double getMeatSubstitutionProportion(){
+ return meatSubstitutionProportion;
+ }
+
public double getSeedWasteRate() {
// just waste rate for animal products
return isAnimalProduct ? ModelConfig.SEED_AND_WASTE_FRACTION/2 : ModelConfig.SEED_AND_WASTE_FRACTION;
--
GitLab