Bug fixes (minimumLandCover FPE).

GAMS/IntExtOpt.gms

@@ -10,7 +10,7 @@
  SET import_crop(all_types)   / monogastrics, ruminants,           wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, sugar, energycrops/;
 
  SET land_cover / cropland, pasture, timberForest, carbonForest, unmanagedForest, otherNatural /;
-* SET wood_producing / timberForest, carbonForest, natural /;
+ SET forest(land_cover) / timberForest, carbonForest, unmanagedForest /;
  ALIAS (land_cover, land_cover_before);
  ALIAS (land_cover, land_cover_after);
 
@@ -120,6 +120,8 @@ $gdxin
 *PARAMETER lcTransitionCost(lc_type_previous, lc_type, location) cost of land transitions including carbon cost and conversion cost;
 *lcTransitionCost(lc_type_previous, lc_type, location) = carbonFluxRate(lc_type_previous, lc_type, location) * carbonPrice + agriExpansionCost(location);
 
+ SCALAR forestExpansionMaxRate / 0.01 /;
+
  VARIABLES
        cropArea(crop, location)             total area for crops - Mha
        fertI(crop, location)              fertilizer intensity for each crop - factor between 0 and 1
@@ -181,6 +183,7 @@ $gdxin
        LAND_COVER_CHANGE_CALC(land_cover, location)
        LAND_COVER_CHANGE_CONSTRAINT(land_cover, location) can't convert more land than was previously available
        LAND_COVER_SELF_CHANGE_CONSTRAINT(land_cover, location) calculate no change area (e.g. cropland to cropland)
+       FOREST_EXPANSION_CONSTRAINT(forest, location)
        WOOD_HARVEST_CALC(location)                          calc wood harvested
        CARBON_FLUX_CALC(location)                           calc carbon flux
        COST_EQ                                        total cost objective function;
@@ -241,19 +244,19 @@ $gdxin
 
  PASTURE_LAND_COVER_CALC(location) .. landCoverArea('pasture', location) =E= cropArea('pasture', location);
 
- MINIMUM_LAND_COVER_CONSTRAINT(land_cover, location) .. landCoverArea(land_cover, location) =G= minimumLandCoverArea(land_cover, location);
-
  LAND_COVER_CHANGE_CALC(land_cover, location) .. landCoverArea(land_cover, location) =E= previousLandCoverArea(land_cover, location) +
                      sum(land_cover_before, landCoverChange(land_cover_before, land_cover, location)) - sum(land_cover_after, landCoverChange(land_cover, land_cover_after, location));
 
-* LAND_COVER_CHANGE_CALC(land_cover, location) .. landCoverArea(land_cover, location) =E= sum(land_cover_before, landCoverChange(land_cover_before, land_cover, location)) -
-*                                                                                         sum(land_cover_after, landCoverChange(land_cover, land_cover_after, location));
-
+ LAND_COVER_CHANGE_CONSTRAINT(land_cover, location) .. sum(land_cover_after, landCoverChange(land_cover, land_cover_after, location)) =L= previousLandCoverArea(land_cover, location);
 
+ LAND_COVER_SELF_CHANGE_CONSTRAINT(land_cover, location) .. landCoverChange(land_cover, land_cover, location) =E= previousLandCoverArea(land_cover, location) -
+                                                            sum(land_cover_after$[not sameAs(land_cover, land_cover_after)], landCoverChange(land_cover, land_cover_after, location));
 
- LAND_COVER_CHANGE_CONSTRAINT(land_cover, location) .. sum(land_cover_after, landCoverChange(land_cover, land_cover_after, location)) =L= previousLandCoverArea(land_cover, location);
+ FOREST_EXPANSION_CONSTRAINT(forest, location) .. landCoverArea(forest, location) - previousLandCoverArea(forest, location) =L= suitableLandArea(location) * forestExpansionMaxRate;
 
- LAND_COVER_SELF_CHANGE_CONSTRAINT(land_cover, location) .. landCoverChange(land_cover, land_cover, location) =E= previousLandCoverArea(land_cover, location) - sum(land_cover_after$[not sameAs(land_cover, land_cover_after)], landCoverChange(land_cover, land_cover_after, location));
+* MINIMUM_LAND_COVER_CONSTRAINT(land_cover, location) .. landCoverArea(land_cover, location) =G= minimumLandCoverArea(land_cover, location);
+* minimum of minimumLandCoverArea and previousLandCoverArea corrects for floating point errors
+  MINIMUM_LAND_COVER_CONSTRAINT(land_cover, location) .. landCoverChange(land_cover, land_cover, location) =G= min(minimumLandCoverArea(land_cover, location), previousLandCoverArea(land_cover, location));
 
  WOOD_HARVEST_CALC(location) .. woodHarvest(location) =E= sum((land_cover_before, land_cover_after), landCoverChange(land_cover_before, land_cover_after, location) * woodYield(land_cover_before, land_cover_after, location));
 

src/ac/ed/lurg/country/CountryAgent.java

@@ -157,7 +157,7 @@ public class CountryAgent extends AbstractCountryAgent {
 		// some hacky code for debug purposes that keeps each gams gdx file
 		try {
 			Files.copy(
-					FileSystems.getDefault().getPath(ModelConfig.TEMP_DIR + File.separator + "_gams_java_gdb1.gdx"),
+					FileSystems.getDefault().getPath(ModelConfig.TEMP_DIR + File.separator + "_gams_java_gdb2.gdx"),
 					FileSystems.getDefault().getPath(ModelConfig.TEMP_DIR + File.separator + country.getName().replaceAll("\\s+","") + currentTimestep.getYear() + ext + ".gdx"),
 					StandardCopyOption.REPLACE_EXISTING);
 		} catch (IOException e) {

src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java

@@ -181,7 +181,7 @@ public class GamsLocationOptimiser {
 				Vector<String> v = new Vector<String>();
 				v.add(lc.getName());
 				v.add(locString);
-				setGamsParamValueTruncate(prevLandCoverP.addRecord(v), landUseItem.getUnprotectedLandCoverArea(lc), 3);
+				setGamsParamValueTruncate(prevLandCoverP.addRecord(v), landUseItem.getUnprotectedLandCoverArea(lc), 5);
 			}
 			
 			/*
@@ -419,7 +419,7 @@ public class GamsLocationOptimiser {
 				Vector<String> v = new Vector<String>();
 				v.add(coverMap.getKey().getName());
 				v.add(locString);
-				setGamsParamValue(minimumLandCover.addRecord(v), coverMap.getValue(), 3);
+				setGamsParamValueTruncate(minimumLandCover.addRecord(v), coverMap.getValue(), 5);
 			}
 		}
 		

src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java

@@ -219,6 +219,11 @@ public class GamsRasterOptimiser {
 
 			}	
 			
+			if (locId == 31) {
+				int foo = 1;
+			}
+				
+			
 			// Allocate land cover change to grid cells
 			for (Entry<LandCoverType, Map<LandCoverType, Double>> fromMap: landCoverChange.getMap().entrySet()) {
 				LandCoverType fromLC = fromMap.getKey();