diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms index a55d1d40e50d35bfb57a5d72450495c34b17bda3..ae12d910cab85da37eae5cddcb387c9ac246fecd 100644 --- a/GAMS/IntExtOpt.gms +++ b/GAMS/IntExtOpt.gms @@ -3,7 +3,7 @@ SET crop(all_types) / wheat, maize, rice, tropicalCereals, oilcrops, soybean, pulses, starchyRoots, pasture /; SET crop_less_pasture(crop) / wheat, maize, rice, tropicalCereals, oilcrops, soybean, pulses, starchyRoots /; SET cereal_crop(crop) / wheat, maize, rice, tropicalCereals /; - SET non_cereal_crop(crop) / oilcrops, soybean, pulses, starchyRoots, pasture /; + SET non_cereal_crop(crop) / oilcrops, soybean, pulses, starchyRoots, pasture /; SET feed_crop(crop) / wheat, maize, oilcrops, soybean, pasture/; SET feed_crop_less_pasture(feed_crop) / wheat, maize, oilcrops, soybean /; SET not_feed_crop(crop) / rice, tropicalCereals, pulses, starchyRoots /; @@ -19,7 +19,7 @@ PARAMETER fertParam(crop, location) yield response to fertilizer parameter; PARAMETER irrigParam(crop, location) yield response to irrigation parameter; PARAMETER demand(all_types) in t; - PARAMETER worldInputEnergy(all_types) average input energy from world exports used to determine if we should import or export; + PARAMETER worldInputEnergy(import_crop) average input energy from world exports used to determine if we should import or export; PARAMETER maxNetImport(import_crop) maximum net import for each crop based on world market; PARAMETER minNetImport(import_crop) minimum net import for each crop based on world market; PARAMETER irrigCost(location) irrigation cost; @@ -40,7 +40,7 @@ $load fertParam, irrigParam, worldInputEnergy, maxNetImport, minNetImport $load meatEfficency, maxLandUseChange, minFeedRate, irrigCost $load cropAdj, fitToPreviousAreas $gdxin - + SCALAR delta use to smooth power function see 7.5 www.gams.com dd docs solversconopt.pdf / 0.00000000001 / SCALAR minDemandPerCereal / 0.1 / @@ -62,7 +62,7 @@ $gdxin soybean 0.5 pulses 0.4 starchyRoots 1.0 - pasture 0.2 / ; + pasture 0.1 / ; SCALAR fertiliserUnitCost / 0.4 / @@ -75,9 +75,10 @@ $gdxin yield(crop, location) yield per area for each crop - t per ha unitEnergy(crop, location) energy per area for each crop - energy net_supply(crop) supply after exports and feed + agriLandExpansion(location) addition agricultural land needed, as it must be positive it deliberately does not account for abandonment energy total input energy - energy; - POSITIVE VARIABLE area, fertI, irrigI, feedAmount; + POSITIVE VARIABLE area, fertI, irrigI, feedAmount, agriLandExpansion; fertI.L(crop, location) = 0.5; irrigI.L(crop, location) = 0.5; @@ -104,6 +105,7 @@ $gdxin PASTURE_IMPORT_CONSTRAINT MIN_FEED_CONSTRAINT constraint on min feed rate NET_SUPPLY_EQ(crop) calc net supply for crops + AGRI_LAND_EXPANSION_CALC(location) calc agriLandExpansion ENERGY_EQ total energy objective function; UNIT_ENERGY_EQ(crop, location) .. unitEnergy(crop, location) =E= baseCost(crop) + @@ -152,14 +154,17 @@ $gdxin MAX_NET_IMPORT_CONSTRAINT(import_crop) .. netImportAmount(import_crop) =L= maxNetImport(import_crop); MIN_NET_IMPORT_CONSTRAINT(import_crop) .. netImportAmount(import_crop) =G= minNetImport(import_crop); - PASTURE_IMPORT_CONSTRAINT .. netImportAmount('pasture') =E= 0; MIN_FEED_CONSTRAINT .. sum(feed_crop_less_pasture, feedDM(feed_crop_less_pasture) * feedAmount(feed_crop_less_pasture)) =G= minFeedRate * (demand('meat') - netImportAmount('meat')); + AGRI_LAND_EXPANSION_CALC(location) .. agriLandExpansion(location) =G= sum(crop, area(crop, location) - previousArea(crop, location)); + +* CROP_TO_PASTURE_CONV_CALC(location) .. cropToPasture(location) =E= sum(crop_less_pasture, area(crop_less_pasture, location)) - area('pasture', location) + ENERGY_EQ .. energy =E= (SUM((crop, location), area(crop, location) * unitEnergy(crop, location)) - + (sum(location, sum(crop, area(crop, location) - previousArea(crop, location))) * landChangeEnergy) - + sum(import_crop, (netImportAmount(import_crop)) * worldInputEnergy(import_crop))) / 1000000; + + (sum(location, agriLandExpansion(location)) * landChangeEnergy) + + sum(import_crop, (netImportAmount(import_crop)) * worldInputEnergy(import_crop))) / 1000000; MODEL LAND_USE /ALL/ ; @@ -169,7 +174,7 @@ $gdxin scalar count / 1 /; scalar newPasture / 0 /; scalar newCropland / 0 /; - scalar pastureAdjRate / 0.5 /; + scalar pastureAdjRate / 0.6 /; * if fitToPreviousAreas is zero then we just optimise if (fitToPreviousAreas eq 0, @@ -194,7 +199,7 @@ $gdxin ); ); - display net_supply.l, demand, feedAmount.l, area.l, yield.l, netImportAmount.l; + display agriLandExpansion.l, previousArea, area.l, net_supply.l, demand, feedAmount.l, yield.l, netImportAmount.l; parameter totalProd(all_types); totalProd(crop) = sum(location, area.l(crop, location) * yield.l(crop, location)); diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java index aeedbbbf9266468a472e8250af323e9b2cc1996e..f94e9f19a0b84095be4e3ab6804a8e753136befe 100644 --- a/src/ac/ed/lurg/ModelMain.java +++ b/src/ac/ed/lurg/ModelMain.java @@ -58,7 +58,10 @@ public class ModelMain { demandManager = new DemandManager(ModelFitType.LOGISTIC, ModelConfig.SSP_SCENARIO); countryAgents = createCountryAgents(); - prevWorldInputCost = new HashMap<CropType, Double>(); // in first timestep we don't have this info, but ok as constrained to import/export specified amount + // in first timestep we don't have this info, but ok as constrained to import/export specified amount + prevWorldInputCost = new HashMap<CropType, Double>(); + for (CropType c : CropType.getImportedTypes()) + prevWorldInputCost.put(c, 0.25); } /* run the model */ @@ -127,7 +130,7 @@ public class ModelMain { } prevWorldInputCost = totalWorldInputCost.divideBy(totalQuantity); - + // output results outputTimestepResults(timestep, globalIntensityRaster, globalCropAreaRaster, globalLocationIdRaster); } @@ -236,7 +239,7 @@ public class ModelMain { // DEBUG code - if (!(country.getCountryName().equals("United States of America") || country.getCountryName().equals("China"))) { + if (!(country.getCountryName().equals("United States of America") )) { //|| country.getCountryName().equals("China") continue; } diff --git a/src/ac/ed/lurg/country/gams/GamsCountryInput.java b/src/ac/ed/lurg/country/gams/GamsCountryInput.java index 7cd568c76b2106594db637330d55df3662cedb34..51e5f44c0d6cb733fc2814649fd32ba5fc35656c 100644 --- a/src/ac/ed/lurg/country/gams/GamsCountryInput.java +++ b/src/ac/ed/lurg/country/gams/GamsCountryInput.java @@ -62,7 +62,7 @@ public class GamsCountryInput { minNetImport = new HashMap<CropType, Double>(); for (Map.Entry<CropType, Double> entry : baseNetImport.entrySet()) { CropType c = entry.getKey(); - double change = allowedImportChange * maxOfProdOrSupply.get(c); + double change = allowedImportChange * maxOfProdOrSupply.get(c); maxNetImport.put(c, entry.getValue() + change); minNetImport.put(c, entry.getValue() - change); } diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java index 500f1f4a42ea4fa4515f9c434aa9e5bc99ffafaa..2699ed0f2d87795196457bb905fd0496a7b3ab74 100644 --- a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java +++ b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java @@ -95,7 +95,7 @@ public class GamsLocationOptimiser { for (CropType cropType : CropType.values()) { double d = areasItem.getCropArea(cropType); - if (DEBUG) LogWriter.println(String.format(" %d %15s,\t %.1f", locationId, cropType.getGamsName(), d)); + if (DEBUG) if (d != 0) LogWriter.println(String.format(" %d %15s,\t %.1f", locationId, cropType.getGamsName(), d)); Vector<String> v = new Vector<String>(); v.add(cropType.getGamsName()); v.add(locString); diff --git a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java index 73f2815240769ed8b90893f8b97cd429ecf99c92..878f10c98a02c11ca06926abb1018d66677c8b1f 100644 --- a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java +++ b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java @@ -199,7 +199,7 @@ public class GamsRasterOptimiser { } - int numCerealCats = 3; + int numCerealCats = 2; int numPastureCats = 2; int thisShouldLookAtCropsOtherThanJustWheat; // need to consider other crops, and perhaps other yieldTypes as well