diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms index 7874b8a2fa104e4de8e0686ce50bb35dbd98d244..7b75fba9f01ac3cb2d13a00c2d5e64dbe8b1ff78 100644 --- a/GAMS/IntExtOpt.gms +++ b/GAMS/IntExtOpt.gms @@ -3,7 +3,7 @@ SET crop_less_pasture(crop) arable crops types includes treenuts but not pasture / cereals, fruits, oilcrops, starchyRoots, treenuts, vegetables, pulses /; SET feed_crop(crop) / cereals, oilcrops /; SET not_feed_crops(crop) / fruits, starchyRoots, treenuts, vegetables, pulses, pasture_or_meat /; - SET location / 1*2 /; + SET location / 1*1 /; PARAMETER previous_area(crop, location) areas for previous timestep in ha; PARAMETER yield00(crop, location) yield in t per ha; @@ -20,13 +20,17 @@ SCALAR tradeBarrier trade barrier which adjust energy cost of imports; SCALAR landChangeEnergy energy required to add ha of agricultural land; SCALAR minFeedRate minimum rate of feed for producing animal products; + SCALAR alpha yield response to fertilizer parameter; + SCALAR beta yield response to irrigation parameter; $gdxin %gdxincname% -$load previous_area, demand, yield00, yield10, yield01, yield11, world_input_energy, maxNetImport, minNetImport, meatEfficency, maxLandUseChange, tradeBarrier, landChangeEnergy, minFeedRate +$load previous_area, demand, yield00, yield10, yield01, yield11, alpha, beta, world_input_energy, maxNetImport, minNetImport, meatEfficency, maxLandUseChange, tradeBarrier, landChangeEnergy, minFeedRate $gdxin DISPLAY yield00, yield10, yield01, yield11 + SCALAR delta use to smooth power function see 7.5 www.gams.com dd docs solversconopt.pdf / 0.000000001 / + PARAMETER feedEnergy(crop) energy from feed in MJ per t / cereals 13.7 oilcrops 18.6 @@ -48,7 +52,7 @@ DISPLAY yield00, yield10, yield01, yield11 unit_energy(crop, location) energy per area for each crop - energy energy total input energy - energy; - POSITIVE VARIABLE area, fertI, irrigI, otherI, feedAmount; + POSITIVE VARIABLE area, fertI, irrigI, otherIntensity, feedAmount; EQUATIONS UNIT_ENERGY_EQ(crop, location) energy per area @@ -59,7 +63,8 @@ DISPLAY yield00, yield10, yield01, yield11 MAX_IRRIG_INTENSITY_CONSTRAINT(crop, location) constraint on maximum irrigation intensity CROPLAND_CHANGE_CONSTRAINT(location) constraint on the rate of land use change PASTURE_CHANGE_CONSTRAINT(location) constraint on the rate of land use change - AGRI_LAND_CHANGE_CONSTRAINT(location) constraint expansion of agricultural land + AGRI_LAND_INCREASE_CONSTRAINT(location) constraint expansion of agricultural land + AGRI_LAND_DECREASE_CONSTRAINT(location) constraint expansion of agricultural land NON_FEED_CROP_CONSTRAINT(not_feed_crops) constraint to set non feed crop feed usage to zero MAX_NET_IMPORT_CONSTRAINT(crop) constraint on max net imports MIN_NET_IMPORT_CONSTRAINT(crop) constraint on min net imports @@ -70,10 +75,10 @@ DISPLAY yield00, yield10, yield01, yield11 sqr(fertI(crop, location)+irrigI(crop, location)+otherIntensity(crop, location)) / 2; YIELD_EQ(crop, location) .. yield(crop, location) =E= ( - yield00(crop, location)*(1-fertI(crop, location))*(1-irrigI(crop, location)) + - yield10(crop, location)*fertI(crop, location)*(1-irrigI(crop, location)) + - yield01(crop, location)*(1-fertI(crop, location))*irrigI(crop, location) + - yield11(crop, location)*fertI(crop, location)*irrigI(crop, location) + yield00(crop, location) + + (yield10(crop, location) - yield00(crop, location)) * ((fertI(crop, location) + delta) ** alpha) + + (yield10(crop, location) - yield00(crop, location)) * ((irrigI(crop, location) + delta) ** beta) + + (yield11(crop, location) + yield00(crop, location) - yield10(crop, location) - yield01(crop, location)) * ((fertI(crop, location) + delta) ** alpha) * ((irrigI(crop, location) + delta) ** beta) ) * otherIntensity(crop, location); CROP_DEMAND_CONSTRAINT(crop_less_pasture) .. sum(location, area(crop_less_pasture, location) * yield(crop_less_pasture, location)) - feedAmount(crop_less_pasture) =G= @@ -90,7 +95,9 @@ DISPLAY yield00, yield10, yield01, yield11 PASTURE_CHANGE_CONSTRAINT(location) .. 1 - area('pasture_or_meat', location)/previous_area('pasture_or_meat', location) =L= maxLandUseChange; - AGRI_LAND_CHANGE_CONSTRAINT(location) .. abs(sum(crop, area(crop, location)) / sum(crop, previous_area(crop, location)) - 1) =L= maxLandUseChange; + AGRI_LAND_INCREASE_CONSTRAINT(location) .. sum(crop, area(crop, location)) / sum(crop, previous_area(crop, location)) =L= 1 + maxLandUseChange; + + AGRI_LAND_DECREASE_CONSTRAINT(location) .. sum(crop, area(crop, location)) / sum(crop, previous_area(crop, location)) =G= 1 - maxLandUseChange; NON_FEED_CROP_CONSTRAINT(not_feed_crops) .. feedAmount(not_feed_crops) =E= 0; @@ -105,7 +112,7 @@ DISPLAY yield00, yield10, yield01, yield11 + sum(crop, (netImportAmount(crop) * tradeBarrier) * world_input_energy(crop)); MODEL LAND_USE /ALL/ ; - SOLVE LAND_USE USING DNLP MINIMIZING energy; + SOLVE LAND_USE USING NLP MINIMIZING energy; DISPLAY yield00, yield10, yield01, yield11 diff --git a/src/ac/ed/lurg/country/CountryAgent.java b/src/ac/ed/lurg/country/CountryAgent.java index 013ba268187dc49ae7474cac4d2f4afb044c26e6..ebc1d8e808a02e46140bf30002bfdad1443ace27 100644 --- a/src/ac/ed/lurg/country/CountryAgent.java +++ b/src/ac/ed/lurg/country/CountryAgent.java @@ -69,10 +69,36 @@ public class CountryAgent implements GamsInputData { } @Override - public Map<CropType, Double> getRefYield() { + public Map<CropType, Double> getYield00() { return currentRefYield; } + @Override + public Map<CropType, Double> getYield10() { + return currentRefYield; + } + + @Override + public Map<CropType, Double> getYield01() { + return currentRefYield; + } + + @Override + public Map<CropType, Double> getYield11() { + return currentRefYield; + } + + @Override + public double getAlpha() { + return 0.4; + } + + @Override + public double getBeta() { + return 0.5; + } + + @Override public Map<CropType, Double> getPreviousCropArea() { Map<CropType, Double> previousCropAreas = new HashMap<CropType, Double>(); diff --git a/src/ac/ed/lurg/country/gams/GamsInputData.java b/src/ac/ed/lurg/country/gams/GamsInputData.java index 3f98c969c04b64b15fa0e2e0b48dd8ebff55303b..bd87cc7dc33a6d3336f69b3bb321c19c4451f5c3 100644 --- a/src/ac/ed/lurg/country/gams/GamsInputData.java +++ b/src/ac/ed/lurg/country/gams/GamsInputData.java @@ -7,7 +7,10 @@ import ac.ed.lurg.types.CropType; public interface GamsInputData { Map<CropType, Double> getProjectedDemand(); - Map<CropType, Double> getRefYield(); + Map<CropType, Double> getYield00(); + Map<CropType, Double> getYield10(); + Map<CropType, Double> getYield01(); + Map<CropType, Double> getYield11(); Map<CropType, Double> getPreviousCropArea(); Map<CropType, Double> getWorldInputEnergy(); Map<CropType, Double> getMaxNetImport(); @@ -18,4 +21,7 @@ public interface GamsInputData { double getTradeBarrier(); double getLandChangeEnergy(); double getMinFeedRate(); + double getAlpha(); + double getBeta(); + } diff --git a/src/ac/ed/lurg/country/gams/GamsLandUseOptimiser.java b/src/ac/ed/lurg/country/gams/GamsLandUseOptimiser.java index bd83bc770b227a22756247048d99baefadeb683a..db8587c9c45ca13d0637b9a4f76b40319b1a91c1 100644 --- a/src/ac/ed/lurg/country/gams/GamsLandUseOptimiser.java +++ b/src/ac/ed/lurg/country/gams/GamsLandUseOptimiser.java @@ -39,10 +39,12 @@ public class GamsLandUseOptimiser { GAMSWorkspace ws = new GAMSWorkspace(wsInfo); GAMSDatabase db = ws.addDatabase(); + int numLocations = 1; + LogWriter.println("\nPrevious areas"); GAMSParameter parm = db.addParameter("previous_area", 2, "the previous area for each land use"); - addLocationMapParm(parm, inputData.getPreviousCropArea(), "1"); - addLocationMapParm(parm, inputData.getPreviousCropArea(), "2"); + for (int i= 1; i<=numLocations; i++) + addLocationMapParm(parm, inputData.getPreviousCropArea(), i); LogWriter.println("\nDemand"); parm = db.addParameter("demand", 1, "demand for crop"); @@ -50,23 +52,23 @@ public class GamsLandUseOptimiser { LogWriter.println("\nYield"); parm = db.addParameter("yield00", 2, "ref yield t per ha"); - addLocationMapParm(parm, inputData.getRefYield(), "1"); - addLocationMapParm(parm, inputData.getRefYield(), "2"); + for (int i= 1; i<=numLocations; i++) + addLocationMapParm(parm, inputData.getYield00(), i); LogWriter.println("\nYield"); parm = db.addParameter("yield10", 2, "ref yield t per ha"); - addLocationMapParm(parm, inputData.getRefYield(), "1"); - addLocationMapParm(parm, inputData.getRefYield(), "2"); + for (int i= 1; i<=numLocations; i++) + addLocationMapParm(parm, inputData.getYield10(), i); LogWriter.println("\nYield"); parm = db.addParameter("yield01", 2, "ref yield t per ha"); - addLocationMapParm(parm, inputData.getRefYield(), "1"); - addLocationMapParm(parm, inputData.getRefYield(), "2"); + for (int i= 1; i<=numLocations; i++) + addLocationMapParm(parm, inputData.getYield01(), i); LogWriter.println("\nYield"); parm = db.addParameter("yield11", 2, "ref yield t per ha"); - addLocationMapParm(parm, inputData.getRefYield(), "1"); - addLocationMapParm(parm, inputData.getRefYield(), "2"); + for (int i= 1; i<=numLocations; i++) + addLocationMapParm(parm, inputData.getYield11(), i); LogWriter.println("\nWorld input energy"); parm = db.addParameter("world_input_energy", 1, "average input energy from world exports used to determine if we should import or export energy per t"); @@ -85,6 +87,8 @@ public class GamsLandUseOptimiser { addScalar(db.addParameter("tradeBarrier", 0, "trade barrier which adjust energy cost of imports"), inputData.getTradeBarrier()); addScalar(db.addParameter("landChangeEnergy", 0, "energy required to add ha of agricultural land"), inputData.getLandChangeEnergy()); addScalar(db.addParameter("minFeedRate", 0, "minimum rate of feed for producing animal products"), inputData.getMinFeedRate()); + addScalar(db.addParameter("alpha", 0), inputData.getAlpha()); + addScalar(db.addParameter("beta", 0), inputData.getBeta()); GAMSJob gamsJob = ws.addJobFromFile(ModelConfig.GAMS_MODEL); GAMSOptions opt = ws.addOptions(); @@ -139,13 +143,13 @@ public class GamsLandUseOptimiser { } } - private void addLocationMapParm(GAMSParameter parm, Map<CropType, Double> itemMap, String locationId) { + private void addLocationMapParm(GAMSParameter parm, Map<CropType, Double> itemMap, int locationId) { for (Map.Entry<CropType, Double> entry : itemMap.entrySet()) { LogWriter.println(String.format(" %15s,\t %.1f", entry.getKey().getGamsName(), entry.getValue())); Vector<String> v = new Vector<String>(); v.add(entry.getKey().getGamsName()); - v.add(locationId); - if (locationId.equals("2") & entry.getKey().getGamsName().equals("cereals")) + v.add(Integer.toString(locationId)); + if (locationId== 2 & entry.getKey().getGamsName().equals("cereals")) parm.addRecord(v).setValue(0.99); else parm.addRecord(v).setValue(entry.getValue()); diff --git a/src/ac/ed/lurg/country/gams/GamsTest.java b/src/ac/ed/lurg/country/gams/GamsTest.java index 94b5901a7c453fcb289afd046932a1d3e49e4935..7a2ed0620e05946d8ce3c98e44a26dc344cd7f3d 100644 --- a/src/ac/ed/lurg/country/gams/GamsTest.java +++ b/src/ac/ed/lurg/country/gams/GamsTest.java @@ -23,28 +23,67 @@ public class GamsTest implements GamsInputData { @Override public Map<CropType, Double> getProjectedDemand() { Map<CropType, Double> dummyMap = new HashMap<CropType, Double>(); - dummyMap.put(CropType.CEREALS, 290.0); - dummyMap.put(CropType.FRUIT, 50.0); + dummyMap.put(CropType.CEREALS, 50.0); + dummyMap.put(CropType.FRUIT, 250.0); dummyMap.put(CropType.OILCROPS, 100.0); - dummyMap.put(CropType.STARCHY_ROOTS, 50.0); - dummyMap.put(CropType.PULSES, 60.0); + dummyMap.put(CropType.STARCHY_ROOTS, 250.0); + dummyMap.put(CropType.PULSES, 260.0); dummyMap.put(CropType.VEGETABLES, 20.0); - dummyMap.put(CropType.MEAT_OR_PASTURE, 80.0); - dummyMap.put(CropType.TREENUTS, 5.0); + dummyMap.put(CropType.MEAT_OR_PASTURE, 800.0); + dummyMap.put(CropType.TREENUTS, 50.0); return dummyMap; } @Override - public Map<CropType, Double> getRefYield() { + public Map<CropType, Double> getYield00() { Map<CropType, Double> dummyMap = new HashMap<CropType, Double>(); dummyMap.put(CropType.CEREALS, 1.0); - dummyMap.put(CropType.FRUIT, 5.0); + dummyMap.put(CropType.FRUIT, 1.0); + dummyMap.put(CropType.OILCROPS, 1.0); + dummyMap.put(CropType.STARCHY_ROOTS, 1.0); + dummyMap.put(CropType.PULSES, 1.0); + dummyMap.put(CropType.VEGETABLES, 1.0); + dummyMap.put(CropType.MEAT_OR_PASTURE, 1.0); + dummyMap.put(CropType.TREENUTS, 1.0); + return dummyMap; + } + + public Map<CropType, Double> getYield10() { + Map<CropType, Double> dummyMap = new HashMap<CropType, Double>(); + dummyMap.put(CropType.CEREALS, 3.0); + dummyMap.put(CropType.FRUIT, 3.0); + dummyMap.put(CropType.OILCROPS, 3.0); + dummyMap.put(CropType.STARCHY_ROOTS, 3.0); + dummyMap.put(CropType.PULSES, 3.0); + dummyMap.put(CropType.VEGETABLES, 3.0); + dummyMap.put(CropType.MEAT_OR_PASTURE, 3.0); + dummyMap.put(CropType.TREENUTS, 3.0); + return dummyMap; + } + + public Map<CropType, Double> getYield01() { + Map<CropType, Double> dummyMap = new HashMap<CropType, Double>(); + dummyMap.put(CropType.CEREALS, 2.0); + dummyMap.put(CropType.FRUIT, 2.0); dummyMap.put(CropType.OILCROPS, 2.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, 8.0); - dummyMap.put(CropType.TREENUTS, 5.0); + dummyMap.put(CropType.STARCHY_ROOTS, 2.0); + dummyMap.put(CropType.PULSES, 2.0); + dummyMap.put(CropType.VEGETABLES, 2.0); + dummyMap.put(CropType.MEAT_OR_PASTURE, 2.0); + dummyMap.put(CropType.TREENUTS, 2.0); + return dummyMap; + } + + public Map<CropType, Double> getYield11() { + Map<CropType, Double> dummyMap = new HashMap<CropType, Double>(); + dummyMap.put(CropType.CEREALS, 4.0); + dummyMap.put(CropType.FRUIT, 4.0); + dummyMap.put(CropType.OILCROPS, 4.0); + dummyMap.put(CropType.STARCHY_ROOTS, 4.0); + dummyMap.put(CropType.PULSES, 4.0); + dummyMap.put(CropType.VEGETABLES, 4.0); + dummyMap.put(CropType.MEAT_OR_PASTURE, 4.0); + dummyMap.put(CropType.TREENUTS, 4.0); return dummyMap; } @@ -133,4 +172,14 @@ public class GamsTest implements GamsInputData { public double getMinFeedRate() { return 0.15; } + + @Override + public double getAlpha() { + return 0.4; + } + + @Override + public double getBeta() { + return 0.5; + } }