Skip to content
Snippets Groups Projects
Commit 59737d6c authored by Peter Alexander's avatar Peter Alexander
Browse files

no message

parent 4c2842a7
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
GAMS/IntExtOpt.gms
+ 15
− 17
View file @ 59737d6c
......@@ -32,7 +32,7 @@
SCALAR fitToPreviousAreas controlling parameter 0 means dont fit and 1 is to calibrate the cropAdj in the calibration year
*$gdxin "/Users/peteralexander/Documents/R_Workspace/temp/GamsTest/_gams_java_189513040.gdx"
*$gdxin "/Users/peteralexander/Documents/R_Workspace/temp/GamsTest/_gams_java_gdb1.gdx"
$gdxin %gdxincname%
$load location, suitableLandArea, previousArea, demand, landChangeEnergy
$load yieldNone, yieldFertOnly, yieldIrrigOnly, yieldBoth
......@@ -75,7 +75,7 @@ $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
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, agriLandExpansion;
......@@ -94,11 +94,11 @@ $gdxin
MAX_FERT_INTENSITY_CONSTRAINT(crop, location) constraint on maximum fertilizer intensity
MAX_IRRIG_INTENSITY_CONSTRAINT(crop, location) constraint on maximum irrigation intensity
TOTAL_LAND_CHANGE_CONSTRAINT(location) constraint on the rate of land use change
CROPLAND_CHANGE_CONSTRAINT(location) constraint on the rate of land use change
PASTURE_MAX_CHANGE_CONSTRAINT(location) constraint on the rate of land use change
PASTURE_MIN_CHANGE_CONSTRAINT(location) constraint on the rate of land use change
AGRI_LAND_INCREASE_CONSTRAINT(location) constraint expansion of agricultural land
AGRI_LAND_DECREASE_CONSTRAINT(location) constraint expansion of agricultural land
* CROPLAND_CHANGE_CONSTRAINT(location) constraint on the rate of land use change
* PASTURE_MAX_CHANGE_CONSTRAINT(location) constraint on the rate of land use change
* PASTURE_MIN_CHANGE_CONSTRAINT(location) constraint on the rate of land use change
* 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_crop) constraint to set non feed crop feed usage to zero
MAX_NET_IMPORT_CONSTRAINT(import_crop) constraint on max net imports
MIN_NET_IMPORT_CONSTRAINT(import_crop) constraint on min net imports
......@@ -142,13 +142,11 @@ $gdxin
TOTAL_LAND_CHANGE_CONSTRAINT(location) .. suitableLandArea(location) =G= sum(crop, area(crop, location));
CROPLAND_CHANGE_CONSTRAINT(location) .. sum(crop_less_pasture, area(crop_less_pasture, location)) =L= (1 + maxLandUseChange) * sum(crop_less_pasture, previousArea(crop_less_pasture, location));
PASTURE_MAX_CHANGE_CONSTRAINT(location) .. area('pasture', location) =L= (1 + maxLandUseChange) * previousArea('pasture', location);
PASTURE_MIN_CHANGE_CONSTRAINT(location) .. area('pasture', location) =G= (1 - maxLandUseChange) * previousArea('pasture', location);
AGRI_LAND_INCREASE_CONSTRAINT(location) .. sum(crop, area(crop, location)) =L= (1 + maxLandUseChange) * sum(crop, previousArea(crop, location));
AGRI_LAND_DECREASE_CONSTRAINT(location) .. sum(crop, area(crop, location)) =G= (1 - maxLandUseChange) * sum(crop, previousArea(crop, location));
* CROPLAND_CHANGE_CONSTRAINT(location) .. sum(crop_less_pasture, area(crop_less_pasture, location)) =L= (1 + maxLandUseChange) * sum(crop_less_pasture, previousArea(crop_less_pasture, location));
* PASTURE_MAX_CHANGE_CONSTRAINT(location) .. area('pasture', location) =L= (1 + maxLandUseChange) * previousArea('pasture', location);
* PASTURE_MIN_CHANGE_CONSTRAINT(location) .. area('pasture', location) =G= (1 - maxLandUseChange) * previousArea('pasture', location);
* AGRI_LAND_INCREASE_CONSTRAINT(location) .. sum(crop, area(crop, location)) =L= (1 + maxLandUseChange) * sum(crop, previousArea(crop, location));
* AGRI_LAND_DECREASE_CONSTRAINT(location) .. sum(crop, area(crop, location)) =G= (1 - maxLandUseChange) * sum(crop, previousArea(crop, location));
NON_FEED_CROP_CONSTRAINT(not_feed_crop) .. feedAmount(not_feed_crop) =E= 0;
......@@ -185,7 +183,7 @@ $gdxin
basePasture = sum(location, previousArea("pasture", location));
baseCropland = sum(location, sum(crop_less_pasture, previousArea(crop_less_pasture, location)));
while((count le 10 and not (abs(basePasture-newPasture) le 0.2 and abs(baseCropland-newCropland) le 0.2)),
while((count le 10 and not (abs(basePasture-newPasture) le 0.1 and abs(baseCropland-newCropland) le 0.1)),
SOLVE LAND_USE USING NLP MINIMIZING energy;
newPasture = sum(location, area.l("pasture", location));
......@@ -213,11 +211,11 @@ $gdxin
totalImportEnergy(import_crop) = netImportAmount.l(import_crop) * worldInputEnergy(import_crop);
feedEnergy(feed_crop)$(netImportAmount.l(feed_crop) gt 0) = (totalProdCost(feed_crop) + totalImportEnergy(feed_crop)) * feedAmount.l(feed_crop) / (totalProd(feed_crop) + netImportAmount.l(feed_crop));
feedEnergy(feed_crop)$(netImportAmount.l(feed_crop) le 0) = totalProdCost(feed_crop) * feedAmount.l(feed_crop) / totalProd(feed_crop);
feedEnergy(feed_crop)$(netImportAmount.l(feed_crop) le 0 and totalProd(feed_crop) gt 0) = totalProdCost(feed_crop) * feedAmount.l(feed_crop) / totalProd(feed_crop);
totalProdCost('meat') = sum(feed_crop, feedEnergy(feed_crop))
display totalProd, totalImportEnergy, feedEnergy, totalProdCost
display totalProdCost, totalProd, totalImportEnergy, feedEnergy, totalProdCost
Scalar ms 'model status', ss 'solve status';
......
src/ac/ed/lurg/ModelMain.java
+ 16
− 14
View file @ 59737d6c
......@@ -2,6 +2,7 @@ package ac.ed.lurg;
import java.io.File;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.Collection;
import java.util.HashMap;
import java.util.HashSet;
......@@ -92,7 +93,7 @@ public class ModelMain {
LogWriter.println("Country " + ca.getCountry());
Collection<RasterKey> countryKeys = countryToKeysMap.get(ca.getCountry());
YieldRaster countryYieldSurfaces = yieldSurfaces.getSubsetRasterForKeys(countryKeys);
// do the optimization
GamsRasterOutput result = ca.determineProduction(timestep, countryYieldSurfaces, prevWorldInputCost);
......@@ -149,11 +150,7 @@ public class ModelMain {
@Override
public int convertToPixelValue(double value) {
if (value == 4.0)
return 10;
else return 1;
//return (int) (value * 3 + 1);
return (int) (value * 3 + 1);
}
}.writeOutput(ModelConfig.WRITE_JPEG_IMAGES);
......@@ -172,26 +169,31 @@ public class ModelMain {
@Override
public int convertToPixelValue(double value) {
return (int) (value * 20);
return (int) (value * 20) + 1;
}
}.writeOutput(ModelConfig.WRITE_JPEG_IMAGES);
new RasterOutputer<AreasItem>(cropAreaRaster, "wheatArea" + timestep) {
outputAreas(timestep, cropAreaRaster, CropType.WHEAT);
outputAreas(timestep, cropAreaRaster, CropType.MAIZE);
outputAreas(timestep, cropAreaRaster, CropType.PASTURE);
}
private void outputAreas(int timestep, RasterSet<AreasItem> cropAreaRaster, final CropType crop) {
new RasterOutputer<AreasItem>(cropAreaRaster, crop.getGamsName() + "Area" + timestep) {
@Override
public Double getValue(RasterKey location) {
AreasItem area = results.get(location);
if (area == null)
return null;
return area.getCropArea(CropType.WHEAT);
return area.getCropArea(crop);
}
@Override
public int convertToPixelValue(double value) {
return (int) (value * 100);
return value > 0 ? 10 : 1;
}
}.writeOutput(ModelConfig.WRITE_JPEG_IMAGES);
}
public Map<Country, List<RasterKey>> createCountryToKeysMap() {
......@@ -231,7 +233,7 @@ public class ModelMain {
RasterSet<IrrigationCostItem> allIrrigationCosts = getIrrigationCosts();
Map<Country, Map<CropType, CropUsageData>> cropUsageDataMap = CropUsageData.readCommodityData();
// HashSet<String> countryExclusionList = new HashSet<String>(Arrays.asList("French Polynesia", "Cabo Verde", "Samoa", "Saint Vincent and the Grenadines"));
HashSet<String> countryExclusionList = new HashSet<String>(Arrays.asList("Sierra Leone", "Togo", "United Arab Emirates")); //"French Polynesia", "Cabo Verde", "Samoa", "Saint Vincent and the Grenadines"));
for (Map.Entry<Country, List<RasterKey>> entry : countryToKeysMap.entrySet()) {
Country country = entry.getKey();
......@@ -239,11 +241,11 @@ 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"))) { //|| country.getCountryName().equals("China")
continue;
}
if (demandManager.getPopulation(country, 2010) < 5) { // countryExclusionList.contains(country.getCountryName()
if (demandManager.getPopulation(country, 2010) < 5 || countryExclusionList.contains(country.getCountryName())) {
LogWriter.printlnError("Skipping " + country);
continue;
}
......
src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
+ 31
− 29
View file @ 59737d6c
......@@ -36,9 +36,9 @@ public class GamsLocationOptimiser {
// oilcrops 0.90
// soy 0.89
// potatoes 0.21
private static final boolean DEBUG = true;
private GamsLocationInput inputData;
public GamsLocationOptimiser(GamsLocationInput inputData) {
......@@ -61,14 +61,14 @@ public class GamsLocationOptimiser {
GAMSJob gamsJob = ws.addJobFromFile(ModelConfig.GAMS_MODEL);
GAMSOptions opt = ws.addOptions();
opt.setAllModelTypes("conopt");
opt.defines("gdxincname", inDB.getName());
long startTime = System.currentTimeMillis();
gamsJob.run(opt, inDB);
// cleanup(ws.workingDirectory());
// cleanup(ws.workingDirectory());
LogWriter.println("Took " + (System.currentTimeMillis() - startTime) + " ms to run");
return handleResults(gamsJob.OutDB());
......@@ -87,7 +87,7 @@ public class GamsLocationOptimiser {
GAMSParameter landP = inDB.addParameter("suitableLandArea", 1);
double totalAgriLand = 0;
for (Map.Entry<Integer, ? extends AreasItem> entry : inputData.getPreviousAreas().entrySet()) {
Integer locationId = entry.getKey();
String locString = Integer.toString(locationId);
......@@ -95,7 +95,8 @@ public class GamsLocationOptimiser {
for (CropType cropType : CropType.values()) {
double d = areasItem.getCropArea(cropType);
if (DEBUG) if (d != 0) LogWriter.println(String.format(" %d %15s,\t %.1f", locationId, cropType.getGamsName(), d));
if (DEBUG && 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);
......@@ -164,16 +165,16 @@ public class GamsLocationOptimiser {
addScalar(inDB.addParameter("meatEfficency", 0), countryInput.getMeatEfficiency());
if (DEBUG) LogWriter.println("\nMeatEfficiency: " + countryInput.getMeatEfficiency());
addScalar(inDB.addParameter("maxLandUseChange", 0), countryInput.getMaxLandUseChange());
if (DEBUG) LogWriter.println("\nMaxLandUseChange: " + countryInput.getMaxLandUseChange());
addScalar(inDB.addParameter("landChangeEnergy", 0), countryInput.getLandChangeEnergy());
if (DEBUG) LogWriter.println("\nLandChangeEnergy: " + countryInput.getLandChangeEnergy());
addScalar(inDB.addParameter("minFeedRate", 0), countryInput.getMinFeedRate());
if (DEBUG) LogWriter.println("\nMinFeedRate: " + countryInput.getMinFeedRate());
double fitToPreviousAreas = countryInput.isCalibrateToObserved() ? 1 : 0;
if (DEBUG) LogWriter.println("\nfitToPreviousAreas: " + fitToPreviousAreas);
addScalar(inDB.addParameter("fitToPreviousAreas", 0), fitToPreviousAreas);
......@@ -191,7 +192,7 @@ public class GamsLocationOptimiser {
GAMSVariable varAreas = outDB.getVariable("area");
GAMSVariable varFertIntensities = outDB.getVariable("fertI");
GAMSVariable varIrrigIntensities = outDB.getVariable("irrigI");
// GAMSVariable varOtherIntensities = outDB.getVariable("otherIntensity");
// GAMSVariable varOtherIntensities = outDB.getVariable("otherIntensity");
GAMSVariable varFeedAmount = outDB.getVariable("feedAmount");
GAMSVariable varNetImports = outDB.getVariable("netImportAmount");
GAMSVariable varYields = outDB.getVariable("yield");
......@@ -199,7 +200,7 @@ public class GamsLocationOptimiser {
GAMSParameter parmCropAdj = outDB.getParameter("cropAdj");
GAMSParameter parmProd = outDB.getParameter("totalProd");
GAMSParameter parmProdCost = outDB.getParameter("totalProdCost");
double totalArea = 0;
double area, fertIntensity, irrigIntensity, otherIntensity = Double.NaN, feedAmount, netImport, yield, unitEnergy, cropAdj, prod, prodCost;
......@@ -215,7 +216,7 @@ public class GamsLocationOptimiser {
area = rec.getLevel();
fertIntensity = varFertIntensities.findRecord(itemName, locationName).getLevel();
irrigIntensity = varIrrigIntensities.findRecord(itemName, locationName).getLevel();
// otherIntensity = varOtherIntensities.findRecord(itemName, locationName).getLevel();
// otherIntensity = varOtherIntensities.findRecord(itemName, locationName).getLevel();
yield = varYields.findRecord(itemName, locationName).getLevel();
unitEnergy = varUnitEnergies.findRecord(itemName, locationName).getLevel();
......@@ -228,34 +229,35 @@ public class GamsLocationOptimiser {
cropAdj = getParmValue(parmCropAdj, itemName);
prod = getParmValue(parmProd, itemName);
prodCost = getParmValue(parmProdCost, itemName);
cropUsageData.put(cropType, new CropUsageData(feedAmount, netImport, prod, prodCost));
cropAdjs.put(cropType, cropAdj);
if (DEBUG) LogWriter.println(String.format("\n%s:\tfeedAmount= %.1f,\tnetImports= %.3f,\tprod= %.3f,\tprodCost= %.3f,\tcropAdj= %.3f", itemName, feedAmount, netImport, prod, prodCost, cropAdj));
}
if (area > 0) {
if (area > 0) {
if (DEBUG) LogWriter.println(String.format("\t location %s, %s:\tarea= %.1f,\tfert= %.3f,\tirrg= %.3f,\tintensity= %.3f",
locationName, itemName, area, fertIntensity, irrigIntensity, otherIntensity));
locationName, itemName, area, fertIntensity, irrigIntensity, otherIntensity));
IntensitiesItem intensityItem = intensities.get(locId);
if (intensityItem == null) {
intensityItem = new IntensitiesItem();
intensities.put(locId, intensityItem);
}
intensityItem.setIntensity(cropType, new Intensity(fertIntensity, irrigIntensity, otherIntensity, yield, unitEnergy));
}
AreasItem areaItem = cropAreas.get(locId);
if (areaItem == null) {
areaItem = new AreasItem();
cropAreas.put(locId, areaItem);
}
areaItem.setCropArea(cropType, area);
totalArea += area;
AreasItem areaItem = cropAreas.get(locId);
if (areaItem == null) {
areaItem = new AreasItem();
cropAreas.put(locId, areaItem);
}
areaItem.setCropArea(cropType, area);
totalArea += area;
prevCropType = cropType;
}
netImport = varNetImports.findRecord(CropType.MEAT.getGamsName()).getLevel();
prod = parmProd.findRecord(CropType.MEAT.getGamsName()).getValue();
prodCost = parmProdCost.findRecord(CropType.MEAT.getGamsName()).getValue();
......@@ -264,11 +266,11 @@ public class GamsLocationOptimiser {
LogWriter.println(String.format("\n%s:\t\t\t\t\tnetImports= %.3f,\tprod= %.3f,\tprodCost= %.3f", CropType.MEAT.getGamsName(), netImport, prod, prodCost));
LogWriter.println(String.format("\nTotal area= %.1f", totalArea));
}
GamsLocationOutput results = new GamsLocationOutput(modelStatus, intensities, cropAreas, cropUsageData, cropAdjs);
return results;
}
private double getParmValue(GAMSParameter aParm, String itemName) {
try {
double d = aParm.findRecord(itemName).getValue();
......@@ -292,7 +294,7 @@ public class GamsLocationOptimiser {
parm.addRecord(entry.getKey().getGamsName()).setValue(d);
}
}
private void addCommodityMapParm(GAMSParameter parm, Map<CommodityType, Double> itemMap) {
for (Map.Entry<CommodityType, Double> entry : itemMap.entrySet()) {
double d = entry.getValue();
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment