From 0239a3784e8f03e478ff9ba919f57bae6a9bfadb Mon Sep 17 00:00:00 2001
From: Peter Alexander <peter@blackhillock.co.uk>
Date: Wed, 17 Dec 2014 12:02:15 +0000
Subject: [PATCH] no message
---
GAMS/IntExtOpt.gms | 6 +--
src/ac/ed/lurg/ModelMain.java | 51 ++++++++++++-------
src/ac/ed/lurg/country/CountryAgent.java | 3 +-
.../country/gams/GamsLocationOptimiser.java | 9 +++-
src/ac/ed/lurg/landuse/IntensitiesItem.java | 2 +-
src/ac/ed/lurg/landuse/Intensity.java | 14 ++++-
src/ac/ed/lurg/types/CropToDoubleMap.java | 35 +++++++++++++
7 files changed, 94 insertions(+), 26 deletions(-)
create mode 100644 src/ac/ed/lurg/types/CropToDoubleMap.java
diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms
index 013b3b9c..da6e1984 100644
--- a/GAMS/IntExtOpt.gms
+++ b/GAMS/IntExtOpt.gms
@@ -50,7 +50,7 @@ $gdxin
feedAmount(crop) amount of feed for each crop - t
netImportAmount(crop) export of crop - t
yield(crop, location) yield per area for each crop - t per ha
- unit_energy(crop, location) energy per area for each crop - energy
+ unitEnergy(crop, location) energy per area for each crop - energy
energy total input energy - energy;
POSITIVE VARIABLE area, fertI, irrigI, otherIntensity, feedAmount;
@@ -73,7 +73,7 @@ $gdxin
MIN_FEED_CONSTRAINT constraint on min feed rate
ENERGY_EQ total energy objective function;
- UNIT_ENERGY_EQ(crop, location) .. unit_energy(crop, location) =E= power((1+ fertI(crop, location)+irrigI(crop, location)+otherIntensity(crop, location)), 2);
+ UNIT_ENERGY_EQ(crop, location) .. unitEnergy(crop, location) =E= power((1+ fertI(crop, location)+irrigI(crop, location)+otherIntensity(crop, location)), 2);
YIELD_EQ(crop, location) .. yield(crop, location) =E= (
yieldNone(crop, location) +
@@ -118,7 +118,7 @@ $gdxin
MIN_FEED_CONSTRAINT .. sum(feed_crop, feedDM(feed_crop) * feedAmount(feed_crop)) =G= minFeedRate * (demand('pasture_or_meat') - netImportAmount('pasture_or_meat'));
- ENERGY_EQ .. energy =E= (SUM((crop, location), area(crop, location) * unit_energy(crop, location))
+ ENERGY_EQ .. energy =E= (SUM((crop, location), area(crop, location) * unitEnergy(crop, location))
+ sum(location, sum(crop, area(crop, location)) - sum(crop, previous_area(crop, location)) * landChangeEnergy)
+ sum(crop, (netImportAmount(crop) * tradeBarrier) * world_input_energy(crop))) / 1000000;
diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java
index bdaa29e8..ae08cb72 100644
--- a/src/ac/ed/lurg/ModelMain.java
+++ b/src/ac/ed/lurg/ModelMain.java
@@ -12,10 +12,15 @@ import ac.ed.lurg.country.Country;
import ac.ed.lurg.country.CountryAgent;
import ac.ed.lurg.country.CountryBoundaryItem;
import ac.ed.lurg.country.CountryBoundaryReader;
+import ac.ed.lurg.country.gams.GamsRasterOutput;
import ac.ed.lurg.demand.DemandManager;
+import ac.ed.lurg.landuse.AreasItem;
import ac.ed.lurg.landuse.CommodityData;
+import ac.ed.lurg.landuse.IntensitiesItem;
+import ac.ed.lurg.landuse.Intensity;
import ac.ed.lurg.landuse.LandCoverItem;
import ac.ed.lurg.landuse.LandCoverReader;
+import ac.ed.lurg.types.CropToDoubleMap;
import ac.ed.lurg.types.CropType;
import ac.ed.lurg.types.LandDataType;
import ac.ed.lurg.types.ModelFitType;
@@ -31,6 +36,8 @@ public class ModelMain {
private Collection<CountryAgent> countryAgents;
private Map<Country, List<RasterKey>> countryToKeysMap;
private DemandManager demandManager;
+ private Map<CropType, Double> prevWorldInputEnergy;
+
public static void main(String[] args) {
ModelMain theModel = new ModelMain();
@@ -42,8 +49,9 @@ public class ModelMain {
private void setup() {
countryToKeysMap = createCountryToKeysMap();
demandManager = new DemandManager(ModelFitType.LOGISTIC, ModelConfig.SSP_SCENARIO);
-
countryAgents = createCountryAgents();
+
+ prevWorldInputEnergy = new HashMap<CropType, Double>(); // in first timestep we don't have this info, but ok as constrained to import/export specified amount
}
/* run the model */
@@ -55,29 +63,36 @@ public class ModelMain {
private void doTimestep(int timestep) {
YieldRaster yieldSurfaces = getYieldSurfaces(timestep);
LogWriter.println("Timestep " + timestep);
+
+ CropToDoubleMap totalArea = new CropToDoubleMap();
+ CropToDoubleMap totalWorldInputEnergy = new CropToDoubleMap();
for (CountryAgent ca : countryAgents) {
LogWriter.println("Country " + ca.getCountry());
Collection<RasterKey> countryKeys = countryToKeysMap.get(ca.getCountry());
YieldRaster countryYieldSurfaces = yieldSurfaces.getSubsetRasterForKeys(countryKeys);
- ca.determineProduction(timestep, countryYieldSurfaces, getWorldInputEnergy());
+ GamsRasterOutput result = ca.determineProduction(timestep, countryYieldSurfaces, prevWorldInputEnergy);
+
+ RasterSet<AreasItem> cropAreaR = result.getCropAreaRaster();
+ RasterSet<IntensitiesItem> intensityR = result.getIntensityRaster();
+
+ for (Map.Entry<RasterKey, IntensitiesItem> entry : intensityR.entrySet()) {
+ RasterKey key = entry.getKey();
+ AreasItem areas = cropAreaR.get(key);
+ IntensitiesItem intensities= entry.getValue();
+
+ for (CropType crop : CropType.values()) {
+ Double area = areas.getCropArea(crop);
+ Intensity intensity = intensities.getIntensity(crop);
+
+ totalArea.incrementValue(crop, area);
+ totalWorldInputEnergy.incrementValue(crop, area * intensity.getUnitEnergy());
+ }
+ }
}
-
- // examine global trade balance
- }
-
- public Map<CropType, Double> getWorldInputEnergy() {
- Map<CropType, Double> dummyMap = new HashMap<CropType, Double>();
- dummyMap.put(CropType.CEREALS, 9.0);
- dummyMap.put(CropType.FRUIT, 5.0);
- dummyMap.put(CropType.OILCROPS, 10.0);
- dummyMap.put(CropType.STARCHY_ROOTS, 5.0);
- dummyMap.put(CropType.PULSES, 20.0);
- dummyMap.put(CropType.VEGETABLES, 20.0);
- dummyMap.put(CropType.MEAT_OR_PASTURE, 80.0);
- dummyMap.put(CropType.TREENUTS, 5.0);
- return dummyMap;
+
+ prevWorldInputEnergy = totalWorldInputEnergy.divideBy(totalArea);
}
public Map<Country, List<RasterKey>> createCountryToKeysMap() {
@@ -145,4 +160,4 @@ public class ModelMain {
return yieldSurfaces;
}
-}
+}
\ No newline at end of file
diff --git a/src/ac/ed/lurg/country/CountryAgent.java b/src/ac/ed/lurg/country/CountryAgent.java
index 87cdf5a5..f2c95fa5 100644
--- a/src/ac/ed/lurg/country/CountryAgent.java
+++ b/src/ac/ed/lurg/country/CountryAgent.java
@@ -67,7 +67,7 @@ public class CountryAgent {
return country;
}
- public void determineProduction(int timestep, YieldRaster countryYieldSurfaces, Map<CropType, Double> worldInputEnergy) {
+ public GamsRasterOutput determineProduction(int timestep, YieldRaster countryYieldSurfaces, Map<CropType, Double> worldInputEnergy) {
currentTimestep = timestep;
this.countryYieldSurfaces = countryYieldSurfaces;
@@ -82,6 +82,7 @@ public class CountryAgent {
GamsRasterOutput result = opti.run();
resultsTimeseries.put(timestep, result);
+ return result;
}
private Map<CropType, Double> getProjectedDemand() {
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
index 1c5a2cad..71cf0eb9 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
@@ -135,9 +135,11 @@ public class GamsLocationOptimiser {
GAMSVariable varOtherIntensities = outDB.getVariable("otherIntensity");
GAMSVariable varFeedAmount = outDB.getVariable("feedAmount");
GAMSVariable varNetImports = outDB.getVariable("netImportAmount");
+ GAMSVariable varYields = outDB.getVariable("yield");
+ GAMSVariable varUnitEnergies = outDB.getVariable("unitEnergy");
double totalArea = 0;
- double area, fertIntensity, irrigIntensity, otherIntensity, feedAmount, netImport;
+ double area, fertIntensity, irrigIntensity, otherIntensity, feedAmount, netImport, yield, unitEnergy;
Map<Integer, IntensitiesItem> intensities = new HashMap<Integer, IntensitiesItem>();
Map<Integer, AreasItem> cropAreas = new HashMap<Integer, AreasItem>();
@@ -151,6 +153,9 @@ public class GamsLocationOptimiser {
fertIntensity = varFertIntensities.findRecord(itemName, locationName).getLevel();
irrigIntensity = varIrrigIntensities.findRecord(itemName, locationName).getLevel();
otherIntensity = varOtherIntensities.findRecord(itemName, locationName).getLevel();
+ yield = varYields.findRecord(itemName, locationName).getLevel();
+ unitEnergy = varUnitEnergies.findRecord(itemName, locationName).getLevel();
+
int locId = Integer.parseInt(locationName);
CropType cropType = CropType.getForGamsName(itemName);
@@ -171,7 +176,7 @@ public class GamsLocationOptimiser {
intensityItem = new IntensitiesItem();
intensities.put(locId, intensityItem);
}
- intensityItem.setIntensity(cropType, new Intensity(fertIntensity, irrigIntensity, otherIntensity));
+ intensityItem.setIntensity(cropType, new Intensity(fertIntensity, irrigIntensity, otherIntensity, yield, unitEnergy));
AreasItem areaItem = cropAreas.get(locId);
if (areaItem == null) {
diff --git a/src/ac/ed/lurg/landuse/IntensitiesItem.java b/src/ac/ed/lurg/landuse/IntensitiesItem.java
index c7006850..90c21543 100644
--- a/src/ac/ed/lurg/landuse/IntensitiesItem.java
+++ b/src/ac/ed/lurg/landuse/IntensitiesItem.java
@@ -10,7 +10,7 @@ public class IntensitiesItem implements RasterItem {
Map<CropType, Intensity> intensityMap = new HashMap<CropType, Intensity>();
- public Intensity getLandUses(CropType crop) {
+ public Intensity getIntensity(CropType crop) {
return intensityMap.get(crop);
}
diff --git a/src/ac/ed/lurg/landuse/Intensity.java b/src/ac/ed/lurg/landuse/Intensity.java
index daa20a8a..6844daaf 100644
--- a/src/ac/ed/lurg/landuse/Intensity.java
+++ b/src/ac/ed/lurg/landuse/Intensity.java
@@ -4,12 +4,16 @@ public class Intensity {
double fertiliserIntensity;
double irrigationIntensity;
double otherIntensity;
+ double yield;
+ double unitEnergy;
- public Intensity(double fertiliserIntensity, double irrigationIntensity, double otherIntensity) {
+ public Intensity(double fertiliserIntensity, double irrigationIntensity, double otherIntensity, double yield, double unitEnergy) {
super();
this.fertiliserIntensity = fertiliserIntensity;
this.irrigationIntensity = irrigationIntensity;
this.otherIntensity = otherIntensity;
+ this.yield = yield;
+ this.unitEnergy = unitEnergy;
}
public double getFertiliserIntensity() {
@@ -23,4 +27,12 @@ public class Intensity {
public double getOtherIntensity() {
return otherIntensity;
}
+
+ public double getYield() {
+ return yield;
+ }
+
+ public double getUnitEnergy() {
+ return unitEnergy;
+ }
}
diff --git a/src/ac/ed/lurg/types/CropToDoubleMap.java b/src/ac/ed/lurg/types/CropToDoubleMap.java
new file mode 100644
index 00000000..857d57e9
--- /dev/null
+++ b/src/ac/ed/lurg/types/CropToDoubleMap.java
@@ -0,0 +1,35 @@
+package ac.ed.lurg.types;
+
+import java.util.HashMap;
+import java.util.Map;
+
+public class CropToDoubleMap extends HashMap<CropType, Double> {
+
+ private static final long serialVersionUID = -4347578627666847612L;
+
+ public void incrementValue(CropType crop, double inc) {
+ if (containsKey(crop)) {
+ double currentValue = get(crop);
+ put(crop, currentValue + inc);
+ }
+ else {
+ put(crop, inc);
+ }
+ }
+
+ /** parameter map must have all the keys in instance*/
+ public CropToDoubleMap divideBy(HashMap<CropType, Double> valuesToDivideBy) {
+ CropToDoubleMap result = new CropToDoubleMap();
+
+ for (Map.Entry<CropType, Double> entry : entrySet()) {
+ CropType crop = entry.getKey();
+ Double divideBy = valuesToDivideBy.get(crop);
+ if (divideBy == null)
+ throw new RuntimeException("CropToDoubleMap.divideBy parameter has no value for " + crop);
+
+ result.put(crop, entry.getValue()/divideBy);
+ }
+ return result;
+ }
+
+}
--
GitLab