diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms
index 9da402d1201d464a84e4644c3b21d8fd6c0b4137..44d288cc8705552d09108d8985ea14fe78413a32 100644
--- a/GAMS/IntExtOpt.gms
+++ b/GAMS/IntExtOpt.gms
@@ -8,6 +8,7 @@
SET feed_crop(crop) / wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, pasture/;
SET feed_crop_less_pasture(feed_crop) / wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg /;
SET import_crop(all_types) / monogastrics, ruminants, wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, sugar, energycrops/;
+ SET market_crop(crop) / wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, sugar, energycrops/;
SET land_cover / cropland, pasture, timberForest, carbonForest, natural /;
SET forested(land_cover) / timberForest, carbonForest, natural /;
@@ -36,8 +37,8 @@
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 exportPrices(import_crop) prices for exports;
- PARAMETER importPrices(import_crop) prices for imports;
+ PARAMETER exportPrices(all_types) prices for exports;
+ PARAMETER importPrices(all_types) prices for imports;
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 in cost per 1000 Mlitre or Mha for each litre per m2;
@@ -73,6 +74,9 @@
SCALAR woodMaxNetImport;
SCALAR woodMinNetImport;
SCALAR managedForestAreaMaxChangeRate;
+ SCALAR forestEstablishmentCost;
+ SCALAR naturalAreaValue;
+ SCALAR managedWoodHarvest;
*$gdxin "C:\Users\Bart\Documents\PhD\GAMS testing area\_gams_java_gdb1.gdx"
$gdxin %gdxincname%
@@ -83,14 +87,31 @@ $load fertParam, irrigParam, otherIParam, exportPrices, importPrices, maxNetImpo
$load meatEfficency, otherICost, irrigCost, irrigMaxRate, irrigConstraint, fertiliserUnitCost, domesticPriceMarkup, minDemandPerCereal, minDemandPerOilcrop, seedAndWasteRate
$load previousLandCoverArea, carbonFluxRate, woodYieldLuc, carbonPrice, woodPrice
$load conversionCost, woodYieldRota, infrastructureCost, forestRotationPeriod, woodMaxNetImport, woodMinNetImport, managedForestAreaMaxChangeRate, woodDemand
+$load forestEstablishmentCost, naturalAreaValue, managedWoodHarvest
$gdxin
- SCALAR delta "use to smooth power function see 7.5 www.gams.com/dd/docs/solversconopt.pdf" / 0.00000000001 /;
+ SCALAR delta "use to smooth power function see 7.5 www.gams.com dd docs solversconopt.pdf" / 0.00000000001 /;
demand(cereal_crop) = demand('cereals') * minDemandPerCereal(cereal_crop);
demand(oilpulse_crop) = demand('oilcropspulses') * minDemandPerOilcrop(oilpulse_crop);
previousCropArea(crop_less_pasture, location) = previousCropArea(crop_less_pasture, location) * (1.0 - unhandledCropRate);
+
+* PARAMETER marketPrices(market_crop);
+* marketPrices('wheat') = exportPrices('wheat');
+* marketPrices('maize') = exportPrices('maize');
+* marketPrices('rice') = exportPrices('rice');
+* marketPrices('oilcrops') = exportPrices('oilcrops');
+* marketPrices('pulses') = exportPrices('pulses');
+* marketPrices('starchyRoots') = exportPrices('starchyRoots');
+* marketPrices('fruitveg') = exportPrices('fruitveg');
+* marketPrices('sugar') = exportPrices('sugar');
+* marketPrices('energycrops') = exportPrices('energycrops');
+
+ PARAMETER maxExport(import_crop);
+ PARAMETER minExport(import_crop);
+ maxExport(import_crop) = max(-minNetImport(import_crop), 0);
+ minExport(import_crop) = max(-maxNetImport(import_crop), 0);
PARAMETER cropDM(crop) kg DM per kg of feed the conversion from feed to meat is done in the R animal product index. Pasture is in DM terms already
/ wheat 0.87
@@ -120,6 +141,8 @@ $gdxin
otherIntCost(crop) = baseCost(crop)*0.7 + otherICost;
baseCost('pasture') = 0.04;
otherIntCost('pasture') = 0.8 + otherICost;
+
+ SCALAR discountRate / 0.08 /;
VARIABLES
cropArea(crop, location) total area for crops - Mha
@@ -132,68 +155,81 @@ $gdxin
exportAmount(all_types) exports of crops and meat - Mt
yield(crop, location) yield per area for each crop - t per ha
unitCost(crop, location) cost per area for each crop - cost
- net_supply(crop) supply after exports and feed
totalFeedDM total feed dry matter
landCoverArea(land_cover, location) land cover area in Mha
landCoverChange(land_cover_before, land_cover_after, location) land cover change in Mha
- totalConversionCost(location) land cover conversion cost - $1000 per ha
- woodHarvest(location) wood harvested from timber forest rotation - Mt C-eq
+ totalConversionCost(land_cover, location) land cover conversion cost - $1000 per ha
+ newForestPotentialHarvest(location)
+ naturalDeforestedArea(location)
+ woodHarvestRota(location) wood harvested from timber forest rotation - Mt C-eq
+ woodHarvestLuc(location)
woodExported total wood sold - Mt C-eq
woodImported
- woodNetSupply
+ woodSupply
+ woodSold
carbonFlux(location) total carbon flux - Mt C
+ supply(all_types)
+ npv
* A "artificial variable for debugging https://www.gams.com/blog/2017/07/misbehaving-model-infeasible/"
total_cost total cost of domestic supply including net imports;
POSITIVE VARIABLE cropArea, fertI, irrigI, otherIntensity, ruminantFeed, monogastricFeed, importAmount, exportAmount,
- totalFeedDM, landCoverArea, landCoverChange, totalConversionCost, woodHarvest, woodExported, woodImported;
+ totalFeedDM, landCoverArea, landCoverChange, totalConversionCost, woodHarvestRota, woodHarvestLuc, woodExported, woodImported,
+ supply, woodSold, woodSupply, exportAmount, newForestPotentialHarvest, naturalDeforestedArea;
* POSITIVE VARIABLE A;
EQUATIONS
UNIT_COST_EQ(crop, location) cost per area - $1000 per ha or $billion per Mha
YIELD_EQ(crop, location) yield given chosen intensity - tonnes per hectare
- CROP_DEMAND_CONSTRAINT(crop) satisfy demand for individual crops
- TOTAL_CEREAL_DEMAND_CONSTRAINT satisfy demand for combined cereals
- TOTAL_OIL_PULSE_DEMAND_CONSTRAINT satisfy demand for combined oilcrops and pulses
- RUMINANT_DEMAND_CONSTRAINT satisfy demand for ruminant products
- MONOGASTRICS_DEMAND_CONSTRAINT satisfy demand for monogastric products
- TOTAL_NON_PASTURE_FEED_DM_CALC calc total feed dry matter not including pasture
- FEED_MIX_CONSTRAINT(feed_crop_less_pasture) limit amount of feed for each feed crop
MAX_FERT_INTENSITY_CONSTRAINT(crop, location) constraint on maximum fertilizer intensity
MAX_IRRIG_INTENSITY_CONSTRAINT(crop, location) constraint on maximum irrigation intensity
MAX_OTHER_INTENSITY_CONSTRAINT(crop, location)
- SETASIDE_AREA_CALC(location)
+ IRRIGATION_CONSTRAINT(location) constraint on water usage
+
+ SUPPLY_CALC(crop)
+ RUMINANTS_SUPPLY_CALC
+ MONOGASTRICS_SUPPLY_CALC
+ TOTAL_NON_PASTURE_FEED_DM_CALC calc total feed dry matter not including pasture
+ FEED_MIX_CONSTRAINT(feed_crop_less_pasture) limit amount of feed for each feed crop
+ SUPPLY_CONSTRAINT(import_crop)
+
MAX_NET_IMPORT_CONSTRAINT(import_crop) constraint on max net imports
MIN_NET_IMPORT_CONSTRAINT(import_crop) constraint on min net imports
PASTURE_IMPORT_CONSTRAINT constraint to not import pasture
- PASTURE_EXPORT_CONSTRAINT constraint to not export pasture
- IRRIGATION_CONSTRAINT(location) constraint on water usage
- NET_SUPPLY_EQ(crop) calc net supply for crops
+ PASTURE_EXPORT_CONSTRAINT
+
+ SETASIDE_AREA_CALC(location)
CROPLAND_LAND_COVER_CALC(location) cropland area equals sum of all crop areas
PASTURE_LAND_COVER_CALC(location) pasture area (land cover) equals pasture area (land use)
- LAND_COVER_CHANGE_CALC(land_cover, location) calc land cover change
+ LAND_COVER_CHANGE_CALC(land_cover, location) "calc land cover change. landCoverChange(A, A, location) defined as unchanged land cover"
LAND_COVER_CHANGE_CONSTRAINT(land_cover, location) conservation of land area
- CONVERSION_COST(location) cost of land cover conversion
-* MANAGED_FOREST_INCREASE_CONSTRAINT(managed_forest)
-* MANAGED_FOREST_DECREASE_CONSTRAINT(managed_forest)
-
- WOOD_HARVEST_CALC(location) calc wood harvested from timber forest rotation
- WOOD_NET_SUPPLY_CALC
- WOOD_DEMAND_CONSTRAINT
- WOOD_EXPORT_CONSTRAINT constraint on maximum wood export
+ CONVERSION_COST(land_cover, location) cost of land cover conversion
+
+ NEW_FOREST_POTENTIAL_HARVEST(location)
+ WOOD_HARVEST_NATURAL_CALC(location)
+ WOOD_HARVEST_NATURAL_CONSTRAINT(location)
+ WOOD_HARVEST_LUC_CALC(location)
+ WOOD_SUPPLY_CALC
+ WOOD_MARKET_CONSTRAINT
+ WOOD_PRODUCTION_CONSTRAINT
WOOD_MIN_TRADE_CONSTRAINT constraint on minimum wood export
WOOD_MAX_TRADE_CONSTRAINT
+
+ CARBON_FLUX_CALC(location) calc carbon flux
+
+ NPV_CALC Net Present Value function;
- CARBON_FLUX_CALC(location) calc carbon flux
- COST_EQ total cost objective function;
+**************** Crop cost **************************
UNIT_COST_EQ(crop, location) .. unitCost(crop, location) =E= ( baseCost(crop) +
fertiliserUnitCost * fertI(crop, location) +
irrigCost(location) * irrigMaxRate(crop, location) * irrigI(crop, location) +
otherIntCost(crop) * otherIntensity(crop, location)
);
+
+**************** Crop yields *******************
YIELD_EQ(crop, location) .. yield(crop, location) =E= (
yieldNone(crop, location) +
@@ -202,81 +238,94 @@ $gdxin
(yieldBoth(crop, location) + yieldNone(crop, location) - yieldFertOnly(crop, location) - yieldIrrigOnly(crop, location)) *
(1 - exp(-fertParam(crop, location)*fertI(crop, location))) * (1 - exp(-irrigParam(crop, location)*irrigI(crop, location)))
) * (1 - exp(-otherIntensity(crop, location)*otherIParam));
+
+ MAX_FERT_INTENSITY_CONSTRAINT(crop, location) .. fertI(crop, location) =L= 1;
+ MAX_IRRIG_INTENSITY_CONSTRAINT(crop, location) .. irrigI(crop, location) =L= 1;
+ MAX_OTHER_INTENSITY_CONSTRAINT(crop, location) .. otherIntensity(crop, location) =L= 1;
+
+ IRRIGATION_CONSTRAINT(location) .. irrigConstraint(location) * suitableLandArea(location) * (1.0 - unhandledCropRate) =G= sum(crop, irrigMaxRate(crop, location) * irrigI(crop, location) * cropArea(crop, location));
- NET_SUPPLY_EQ(crop) .. net_supply(crop) =E= sum(location, cropArea(crop, location) * yield(crop, location)) * (1 - seedAndWasteRate(crop)) - ruminantFeed(crop) - monogastricFeed(crop) + importAmount(crop) - exportAmount(crop);
-
- CROP_DEMAND_CONSTRAINT(crop) .. net_supply(crop) =G= demand(crop);
-
- TOTAL_CEREAL_DEMAND_CONSTRAINT .. sum(cereal_crop, net_supply(cereal_crop)) =G= demand('cereals');
-
- TOTAL_OIL_PULSE_DEMAND_CONSTRAINT .. sum(oilpulse_crop, net_supply(oilpulse_crop)) =G= demand('oilcropspulses');
-
- RUMINANT_DEMAND_CONSTRAINT .. meatEfficency * sum(feed_crop, ruminantFeed(feed_crop) * cropDM(feed_crop)) * (1 - seedAndWasteRate('ruminants')) =G= (demand('ruminants') - importAmount('ruminants') + exportAmount('ruminants'));
+*************** Commodity supply *************************
- MONOGASTRICS_DEMAND_CONSTRAINT .. meatEfficency * sum(feed_crop_less_pasture, monogastricFeed(feed_crop_less_pasture) * cropDM(feed_crop_less_pasture)) * (1 - seedAndWasteRate('monogastrics')) =G= (demand('monogastrics') - importAmount('monogastrics') + exportAmount('monogastrics'));
+ SUPPLY_CALC(crop) .. supply(crop) =E= sum(location, cropArea(crop, location) * yield(crop, location)) * (1 - seedAndWasteRate(crop)) - ruminantFeed(crop) - monogastricFeed(crop);
+
+ RUMINANTS_SUPPLY_CALC .. supply('ruminants') =E= meatEfficency * sum(feed_crop, ruminantFeed(feed_crop) * cropDM(feed_crop)) * (1 - seedAndWasteRate('ruminants'));
+
+ MONOGASTRICS_SUPPLY_CALC .. supply('monogastrics') =E= meatEfficency * sum(feed_crop_less_pasture, monogastricFeed(feed_crop_less_pasture) * cropDM(feed_crop_less_pasture)) * (1 - seedAndWasteRate('monogastrics'));
TOTAL_NON_PASTURE_FEED_DM_CALC .. totalFeedDM =E= sum(feed_crop_less_pasture, (ruminantFeed(feed_crop_less_pasture) + monogastricFeed(feed_crop_less_pasture)) * cropDM(feed_crop_less_pasture));
FEED_MIX_CONSTRAINT(feed_crop_less_pasture) .. (ruminantFeed(feed_crop_less_pasture) + monogastricFeed(feed_crop_less_pasture)) * cropDM(feed_crop_less_pasture) =L= totalFeedDM * 0.7;
- MAX_FERT_INTENSITY_CONSTRAINT(crop, location) .. fertI(crop, location) =L= 1;
- MAX_IRRIG_INTENSITY_CONSTRAINT(crop, location) .. irrigI(crop, location) =L= 1;
- MAX_OTHER_INTENSITY_CONSTRAINT(crop, location) .. otherIntensity(crop, location) =L= 1;
+ SUPPLY_CONSTRAINT(import_crop) .. supply(import_crop) =E= demand(import_crop) + exportAmount(import_crop) - importAmount(import_crop);
- SETASIDE_AREA_CALC(location) .. cropArea('setaside', location) =E= sum(crop_less_pasture, cropArea(crop_less_pasture, location)) * setAsideRate;
+*************** Exports *******************************
MAX_NET_IMPORT_CONSTRAINT(import_crop) .. importAmount(import_crop) - exportAmount(import_crop) =L= maxNetImport(import_crop);
MIN_NET_IMPORT_CONSTRAINT(import_crop) .. importAmount(import_crop) - exportAmount(import_crop) =G= minNetImport(import_crop);
PASTURE_IMPORT_CONSTRAINT .. importAmount('pasture') =E= 0;
PASTURE_EXPORT_CONSTRAINT .. exportAmount('pasture') =E= 0;
- IRRIGATION_CONSTRAINT(location) .. irrigConstraint(location) * suitableLandArea(location) * (1.0 - unhandledCropRate) =G= sum(crop, irrigMaxRate(crop, location) * irrigI(crop, location) * cropArea(crop, location));
+************** Land Cover *****************************
+ SETASIDE_AREA_CALC(location) .. cropArea('setaside', location) =E= sum(crop_less_pasture, cropArea(crop_less_pasture, location)) * setAsideRate;
+
CROPLAND_LAND_COVER_CALC(location) .. landCoverArea('cropland', location) =E= sum(crop_less_pasture, cropArea(crop_less_pasture, location)) / (1.0 - unhandledCropRate);
-
+
PASTURE_LAND_COVER_CALC(location) .. landCoverArea('pasture', location) =E= cropArea('pasture', 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));
-* landCoverChange(A, A, location) defined as unchanged land cover
+
LAND_COVER_CHANGE_CONSTRAINT(land_cover, location) .. sum(land_cover_after, landCoverChange(land_cover, land_cover_after, location)) =E= previousLandCoverArea(land_cover, location);
- CONVERSION_COST(location) .. totalConversionCost(location) =E= sum((land_cover_before, land_cover_after), landCoverChange(land_cover_before, land_cover_after, location) * conversionCost(land_cover_before, land_cover_after));
+ CONVERSION_COST(land_cover, location) .. totalConversionCost(land_cover, location) =E= sum(land_cover_before, landCoverChange(land_cover_before, land_cover, location) * conversionCost(land_cover_before, land_cover));
-* MANAGED_FOREST_INCREASE_CONSTRAINT(managed_forest) .. sum(location, landCoverArea(managed_forest, location)) =L= sum(location, previousLandCoverArea(managed_forest, location)) * (1 + managedForestAreaMaxChangeRate);
-* MANAGED_FOREST_DECREASE_CONSTRAINT(managed_forest) .. sum(location, landCoverArea(managed_forest, location)) =G= sum(location, previousLandCoverArea(managed_forest, location)) * (1 - managedForestAreaMaxChangeRate);
+************* Forestry ***********************************
-* Wood from timberForest rotation and new timberForest
- WOOD_HARVEST_CALC(location) .. woodHarvest(location) =E= sum(land_cover, landCoverChange(land_cover, 'timberForest', location) * woodYieldRota(land_cover, 'timberForest', location) / forestRotationPeriod) +
- sum((land_cover_before, land_cover_after), landCoverChange(land_cover_before, land_cover_after, location) * woodYieldLuc(land_cover_before, land_cover_after, location));
+ NEW_FOREST_POTENTIAL_HARVEST(location) .. newForestPotentialHarvest(location) =E= sum(land_cover, landCoverChange(land_cover, 'timberForest', location) * woodYieldRota(land_cover, 'timberForest', location));
- WOOD_NET_SUPPLY_CALC .. woodNetSupply =E= sum(location, woodHarvest(location)) - woodExported + woodImported;
+ WOOD_HARVEST_NATURAL_CALC(location) .. woodHarvestNat(location) =E= naturalDeforestedArea(location) * woodYield('natural', location);
- WOOD_DEMAND_CONSTRAINT .. woodNetSupply =G= woodDemand;
+ WOOD_HARVEST_NATURAL_CONSTRAINT(location) .. naturalDeforestedArea(location) =L= landCoverChange('natural', 'natural' location);
+
+ WOOD_HARVEST_LUC_CALC(location) .. woodHarvestLuc(location) =E= sum((forested, land_cover_after), landCoverChange(forested, land_cover_after, location) * woodYieldLuc(forested, land_cover_after, location));
- WOOD_EXPORT_CONSTRAINT .. woodExported =L= sum(location, woodHarvest(location));
+ WOOD_SUPPLY_CALC .. woodSupply =E= managedWoodHarvest + sum(location, woodHarvestNat(location) + woodHarvestLuc(location));
+
+ WOOD_MARKET_CONSTRAINT .. woodSold =E= woodDemand + woodExported - woodImported;
+
+ WOOD_PRODUCTION_CONSTRAINT .. woodSupply =G= woodSold;
WOOD_MIN_TRADE_CONSTRAINT .. woodImported - woodExported =L= woodMaxNetImport;
WOOD_MAX_TRADE_CONSTRAINT .. woodImported - woodExported =G= woodMinNetImport;
+
+*********** Carbon fluxes ***********************************
CARBON_FLUX_CALC(location) .. carbonFlux(location) =E= sum((land_cover_before, land_cover_after), landCoverChange(land_cover_before, land_cover_after, location) * carbonFluxRate(land_cover_before, land_cover_after, location));
+************ Net Present Value ******************************
+
+ NPV_CALC .. npv =E= {
+ [
+ sum(import_crop, exportPrices(import_crop) * supply(import_crop)) +
+ sum(location, landCoverArea('natural', location)) * naturalAreaValue
+ ] * exp(-discountRate) -
+ sum((crop, location), cropArea(crop, location) * unitCost(crop, location) * (1-subsidyRate(crop))) -
+ sum(import_crop, importPrices(import_crop) * importAmount(import_crop))
+ } / (1 - exp(-discountRate) +
+
+ {
+ woodPrice * newForestPotentialHarvest(location) * exp(-discountRate * forestRotationPeriod) -
+ sum(location, landCoverArea('timberForest', location) * forestEstablishmentCost) -
+ woodPrice * woodImported
+ } / (1 - exp(-discountRate * forestRotationPeriod)) -
+
+ sum((land_cover, location), totalConversionCost(land_cover, location));
+
- COST_EQ .. total_cost =E=
- (
- [ SUM((crop, location), cropArea(crop, location) * unitCost(crop, location) * (1-subsidyRate(crop))) ] * domesticPriceMarkup +
-
- SUM(import_crop, importAmount(import_crop) * importPrices(import_crop) - exportAmount(import_crop) * exportPrices(import_crop)) +
-
- SUM(location, totalConversionCost(location)) +
-
- [woodPrice * (woodImported - woodExported)] +
-
- SUM(location, carbonFlux(location)) * carbonPrice
- );
MODEL LAND_USE /ALL/ ;
fertI.L(crop, location) = previousFertIntensity(crop, location);
@@ -289,18 +338,17 @@ $gdxin
monogastricFeed.L(feed_crop) = previousMonogastricFeed(feed_crop);
importAmount.L(all_types) = previousImportAmount(all_types);
exportAmount.L(all_types) = previousExportAmount(all_types);
- totalConversionCost.L(location) = sum(land_cover, landCoverChange.L(land_cover, land_cover, location) * conversionCost(land_cover, land_cover));
- woodHarvest.L(location) = previousLandCoverArea('timberForest', location) * woodYieldRota('timberForest', 'timberForest', location) / forestRotationPeriod;
-
+ totalConversionCost.L(land_cover, location) = landCoverChange.L(land_cover, land_cover, location) * conversionCost(land_cover, land_cover);
* LAND_USE.OptFile = 1;
* display landCoverChange.L
- SOLVE LAND_USE USING NLP MINIMIZING total_cost;
+* SOLVE LAND_USE USING NLP MINIMIZING total_cost;
+ SOLVE LAND_USE USING NLP MAXIMIZING npv;
display previousCropArea, irrigMaxRate, otherIntensity.L, fertI.L, irrigI.L, cropArea.L;
- display net_supply.l, demand, ruminantFeed.l, monogastricFeed.l, importAmount.l, exportAmount.l;
+ display supply.l, demand, ruminantFeed.l, monogastricFeed.l, exportAmount.l;
* Calculate summary information used in Java process
parameter totalProd(all_types);
@@ -331,12 +379,14 @@ $gdxin
totalProdCost('ruminants') = sum(feed_crop, ruminantFeed.l(feed_crop) * feedCostRate(feed_crop));
totalProdCost('monogastrics') = sum(feed_crop, monogastricFeed.l(feed_crop) * feedCostRate(feed_crop));
- netImportCost(import_crop) = importAmount.l(import_crop) * importPrices(import_crop) - exportAmount.l(import_crop) * exportPrices(import_crop);
netImportAmount(import_crop) = importAmount.l(import_crop) - exportAmount.l(import_crop);
+ netImportCost(import_crop) = importAmount.l(import_crop) * importPrices(import_crop) - exportAmount.l(import_crop) * exportPrices(import_crop);
+
netCarbonFlux = SUM(location, carbonFlux.L(location));
netWoodImport = woodImported.L - woodExported.L;
- totalWoodHarvest = sum(location, woodHarvest.L(location));
+ netWoodImport = 0;
+ totalWoodHarvest = sum(location, woodHarvestRota.L(location) + woodHarvestLuc.L(location));
Scalar totalCostsLessLU;
diff --git a/GAMS/IntExtOpt.~gm b/GAMS/IntExtOpt.~gm
new file mode 100644
index 0000000000000000000000000000000000000000..b59dfcf68239b6e5a84400d277207d8c31b21039
--- /dev/null
+++ b/GAMS/IntExtOpt.~gm
@@ -0,0 +1,382 @@
+
+ SET all_types / monogastrics, ruminants, cereals, oilcropspulses, wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, sugar, energycrops, pasture, setaside/;
+
+ SET crop(all_types) / wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, sugar, energycrops, pasture, setaside /;
+ SET crop_less_pasture(crop) / wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, sugar, energycrops, setaside/;
+ SET cereal_crop(crop) / wheat, maize, rice /;
+ SET oilpulse_crop(crop) / oilcrops, pulses /;
+ SET feed_crop(crop) / wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, pasture/;
+ SET feed_crop_less_pasture(feed_crop) / wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg /;
+ SET import_crop(all_types) / monogastrics, ruminants, wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, sugar, energycrops/;
+
+ SET land_cover / cropland, pasture, timberForest, carbonForest, natural/;
+ SET forested(land_cover) / timberForest, carbonForest, natural /;
+ SET managed_forest(land_cover) / timberForest, carbonForest /;
+ SET non_timber_forest(land_cover) / cropland, pasture, carbonForest, natural/;
+ SET exc_natural(land_cover) / cropland, pasture, timberForest/;
+ ALIAS (land_cover, land_cover_before);
+ ALIAS (land_cover, land_cover_after);
+
+ SET location;
+ PARAMETER suitableLandArea(location) areas of land in Mha;
+ PARAMETER previousCropArea(crop, location) areas for previous timestep in Mha;
+ PARAMETER previousFertIntensity(crop, location);
+ PARAMETER previousIrrigIntensity(crop, location);
+ PARAMETER previousOtherIntensity(crop, location);
+ PARAMETER previousRuminantFeed(crop);
+ PARAMETER previousMonogastricFeed(crop);
+ PARAMETER previousImportAmount(all_types);
+ PARAMETER previousExportAmount(all_types);
+ PARAMETER yieldNone(crop, location) yield in t per ha;
+ PARAMETER yieldFertOnly(crop, location) yield in t per ha;
+ PARAMETER yieldIrrigOnly(crop, location) yield in t per ha;
+ PARAMETER yieldBoth(crop, location) yield in t per ha;
+ PARAMETER yieldShock(crop, location) rate of yield shock;
+ 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 exportPrices(import_crop) prices for exports;
+ PARAMETER importPrices(import_crop) prices for imports;
+ 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 in cost per 1000 Mlitre or Mha for each litre per m2;
+ PARAMETER irrigMaxRate(crop, location) max water application rate irrigation in litre per m2;
+ PARAMETER irrigConstraint(location) max water available for irrigation in litre per m2;
+ PARAMETER minDemandPerCereal(cereal_crop) min demand for each cereal crop as factor of all cereals;
+ PARAMETER minDemandPerOilcrop(oilpulse_crop) min demand for oilcrop pulses as factor of total;
+ PARAMETER seedAndWasteRate(all_types) rate of use for seed and waste combined;
+ PARAMETER subsidyRate(crop) rates of subsidy compared to costs;
+
+ SCALAR meatEfficency efficiency of converting feed and pasture into animal products;
+ SCALAR fertiliserUnitCost fert cost at max fert rate;
+ SCALAR otherIParam yield response to other intensity;
+ SCALAR otherICost cost of other intensity;
+ SCALAR unhandledCropRate rate of fruit veg and other crops not modelled;
+ SCALAR setAsideRate rate of set aside and failed crop;
+ SCALAR domesticPriceMarkup factor price increased from cost of production;
+ SCALAR maxLandExpansionRate max rate of country land expansion;
+
+ PARAMETER previousLandCoverArea(land_cover, location) land cover area in Mha;
+
+ PARAMETER carbonFluxRate(land_cover_before, land_cover_after, location) carbon flux - MtC-eq per Mha;
+ PARAMETER woodYieldLuc(land_cover_before, land_cover_after, location) wood yield - MtC-eq per Mha;
+ PARAMETER woodYieldRota(land_cover_before, land_cover_after, location)
+ SCALAR carbonPrice price of carbon - $1000 per tonne;
+ SCALAR woodPrice price of wood and timber - $1000 per tC-eq;
+ PARAMETER timberForestYield(location);
+
+ PARAMETER conversionCost(land_cover, land_cover) cost of converting from one land cover to another;
+ PARAMETER infrastructureCost(location);
+
+ SCALAR netPresentValueFactor;
+ SCALAR forestRotationPeriod;
+ SCALAR timberMaxNetImport;
+ SCALAR timberMinNetImport;
+
+*$gdxin "C:\Users\Bart\Documents\PhD\GAMS testing area\_gams_java_gdb1.gdx"
+$gdxin %gdxincname%
+$load location, suitableLandArea, demand
+$load previousCropArea, previousFertIntensity, previousIrrigIntensity, previousOtherIntensity, previousRuminantFeed, previousMonogastricFeed, previousImportAmount, previousExportAmount
+$load yieldNone, yieldFertOnly, yieldIrrigOnly, yieldBoth, yieldShock
+$load fertParam, irrigParam, otherIParam, exportPrices, importPrices, maxNetImport, minNetImport, unhandledCropRate, setAsideRate, maxLandExpansionRate, subsidyRate
+$load meatEfficency, otherICost, irrigCost, irrigMaxRate, irrigConstraint, fertiliserUnitCost, domesticPriceMarkup, minDemandPerCereal, minDemandPerOilcrop, seedAndWasteRate
+$load previousLandCoverArea, carbonFluxRate, woodYieldLuc, carbonPrice, woodPrice
+$load conversionCost, timberForestYield, woodYieldRota, infrastructureCost, netPresentValueFactor, forestRotationPeriod, timberMaxNetImport, timberMinNetImport
+$gdxin
+
+ SCALAR delta "use to smooth power function see 7.5 www.gams.com/dd/docs/solversconopt.pdf" / 0.00000000001 /;
+
+ demand(cereal_crop) = demand('cereals') * minDemandPerCereal(cereal_crop);
+ demand(oilpulse_crop) = demand('oilcropspulses') * minDemandPerOilcrop(oilpulse_crop);
+
+ previousCropArea(crop_less_pasture, location) = previousCropArea(crop_less_pasture, location) * (1.0 - unhandledCropRate);
+
+ PARAMETER cropDM(crop) kg DM per kg of feed the conversion from feed to meat is done in the R animal product index. Pasture is in DM terms already
+ / wheat 0.87
+ maize 0.86
+ rice 0.89
+ oilcrops 0.96
+ pulses 0.31
+ starchyRoots 0.25
+ fruitveg 0.1
+ sugar 1
+ pasture 1 / ;
+
+ PARAMETER prodCost(crop) cost per ha before intensity values not including fertiliser or irrigation i.e. seed spray and other. In 1000 $ per ha
+ / wheat 0.32
+ maize 0.31
+ rice 0.36
+ oilcrops 0.2
+ pulses 0.31
+ starchyRoots 3.14
+ fruitveg 4.0
+ sugar 3.0
+ energyCrops 0.34 / ;
+
+ PARAMETER baseCost(crop);
+ PARAMETER otherIntCost(crop);
+ baseCost(crop) = prodCost(crop)*0.3;
+ otherIntCost(crop) = baseCost(crop)*0.7 + otherICost;
+ baseCost('pasture') = 0.04;
+ otherIntCost('pasture') = 0.8 + otherICost;
+
+ SCALAR naturalDecreaseMaxRate / 1 /;
+ SCALAR naturalWoodHarvestMaxRate / 0.02 /;
+ woodYieldLuc('carbonForest', 'naturalNew', location) = woodYieldLuc('carbonForest', 'pasture', location);
+ woodYieldLuc('natural', 'naturalNew', location) = woodYieldLuc('natural', 'pasture', location);
+ conversionCost(land_cover, 'naturalNew') = conversionCost(land_cover, 'natural');
+ conversionCost('natural', 'naturalNew') = 0;
+
+ woodYieldLuc('carbonForest', land_cover, location) = 0;
+
+ carbonFluxRate('natural', 'naturalNew', location) = carbonFluxRate('natural', 'pasture', location);
+ carbonFluxRate('carbonForest', 'naturalNew', location) = carbonFluxRate('natural', 'pasture', location);
+ carbonFluxRate('pasture', 'naturalNew', location) = carbonFluxRate('pasture', 'natural', location);
+ carbonFluxRate('cropland', 'naturalNew', location) = carbonFluxRate('cropland', 'natural', location);
+
+ SCALAR forestRotationCost;
+ forestRotationCost = 0.1;
+ SCALAR woodHarvestParam / 1.897 /;
+
+ VARIABLES
+ cropArea(crop, location) total area for crops - Mha
+ fertI(crop, location) fertilizer intensity for each crop - factor between 0 and 1
+ irrigI(crop, location) irrigation intensity for each crop - factor between 0 and 1
+ otherIntensity(crop, location)
+ ruminantFeed(crop) amount of feed for ruminant animals - Mt
+ monogastricFeed(crop) amount of feed for monogatric animals - Mt
+ importAmount(all_types) imports of crops and meat - Mt
+ exportAmount(all_types) exports of crops and meat - Mt
+ yield(crop, location) yield per area for each crop - t per ha
+ unitCost(crop, location) cost per area for each crop - cost
+ net_supply(crop) supply after exports and feed
+ totalFeedDM total feed dry matter
+ landCoverArea(land_cover, location) land cover area in Mha
+ landCoverChange(land_cover_before, land_cover_after, location) land cover change in Mha
+ totalConversionCost(location) land cover conversion cost - $1000 per ha
+ woodHarvestLuc(location) wood harvested due to land cover change - Mt C-eq
+ woodHarvestRota(location) wood harvested from timber forest rotation - Mt C-eq
+ totalWoodHarvest(location) total wood harvested from all activities - Mt C-eq
+ woodExported total wood sold - Mt C-eq
+ carbonFlux(location) total carbon flux - Mt C
+ loggingCost(location) cost of logging forest
+* A "artificial variable for debugging https://www.gams.com/blog/2017/07/misbehaving-model-infeasible/"
+ total_cost total cost of domestic supply including net imports;
+
+ POSITIVE VARIABLE cropArea, fertI, irrigI, otherIntensity, ruminantFeed, monogastricFeed, importAmount, exportAmount,
+ totalFeedDM, landCoverArea, landCoverChange, totalConversionCost, woodHarvestLuc, woodHarvestRota,
+ totalWoodHarvest, woodExported
+
+* POSITIVE VARIABLE A;
+
+ EQUATIONS
+ UNIT_COST_EQ(crop, location) cost per area - $1000 per ha or $billion per Mha
+ YIELD_EQ(crop, location) yield given chosen intensity - tonnes per hectare
+ CROP_DEMAND_CONSTRAINT(crop) satisfy demand for individual crops
+ TOTAL_CEREAL_DEMAND_CONSTRAINT satisfy demand for combined cereals
+ TOTAL_OIL_PULSE_DEMAND_CONSTRAINT satisfy demand for combined oilcrops and pulses
+ RUMINANT_DEMAND_CONSTRAINT satisfy demand for ruminant products
+ MONOGASTRICS_DEMAND_CONSTRAINT satisfy demand for monogastric products
+ TOTAL_NON_PASTURE_FEED_DM_CALC calc total feed dry matter not including pasture
+ FEED_MIX_CONSTRAINT(feed_crop_less_pasture) limit amount of feed for each feed crop
+ MAX_FERT_INTENSITY_CONSTRAINT(crop, location) constraint on maximum fertilizer intensity
+ MAX_IRRIG_INTENSITY_CONSTRAINT(crop, location) constraint on maximum irrigation intensity
+ MAX_OTHER_INTENSITY_CONSTRAINT(crop, location)
+ SETASIDE_AREA_CALC(location)
+ MAX_NET_IMPORT_CONSTRAINT(import_crop) constraint on max net imports
+ MIN_NET_IMPORT_CONSTRAINT(import_crop) constraint on min net imports
+ PASTURE_IMPORT_CONSTRAINT constraint to not import pasture
+ PASTURE_EXPORT_CONSTRAINT constraint to not export pasture
+ IRRIGATION_CONSTRAINT(location) constraint on water usage
+ NET_SUPPLY_EQ(crop) calc net supply for crops
+
+ CROPLAND_LAND_COVER_CALC(location) cropland area equals sum of all crop areas
+ PASTURE_LAND_COVER_CALC(location) pasture area (land cover) equals pasture area (land use)
+ LAND_COVER_CHANGE_CALC(land_cover, location) calc land cover change
+ LAND_COVER_CHANGE_CONSTRAINT(land_cover, location) conservation of land area
+ NATURAL_CONVERSION_CONSTRAINT(exc_natural, location) can't create primary natural area
+ NATURAL_CONTRACTION_CONSTRAINT(location)
+ CONVERSION_COST(location) cost of land cover conversion
+
+ WOOD_HARVEST_LUC_CALC(location) calc wood harvested from LUC
+ WOOD_HARVEST_ROTA_CALC(location) calc wood harvested from timber forest rotation
+ TOTAL_WOOD_HARVEST_CALC(location) total wood harvested
+ NATURAL_WOOD_HARVEST_CONSTRAINT(location)
+ WOOD_MAX_EXPORT_CONSTRAINT constraint on maximum wood export
+ WOOD_MIN_TRADE_CONSTRAINT constraint on minimum wood export
+ WOOD_MAX_TRADE_CONSTRAINT
+ WOOD_ROTA_EXPORT_CONSTRAINT must export all wood produced from timber forest rotation
+
+* WOOD_ROTA_TO_LUC_CONSTRAINT
+ LOGGING_COST_CALC(location)
+ CARBON_FLUX_CALC(location) calc carbon flux
+ COST_EQ total cost objective function;
+
+ UNIT_COST_EQ(crop, location) .. unitCost(crop, location) =E= ( baseCost(crop) +
+ fertiliserUnitCost * fertI(crop, location) +
+ irrigCost(location) * irrigMaxRate(crop, location) * irrigI(crop, location) +
+ otherIntCost(crop) * otherIntensity(crop, location)
+ );
+
+ YIELD_EQ(crop, location) .. yield(crop, location) =E= (
+ yieldNone(crop, location) +
+ (yieldFertOnly(crop, location) - yieldNone(crop, location)) * (1 - exp(-fertParam(crop, location)*fertI(crop, location))) +
+ (yieldIrrigOnly(crop, location) - yieldNone(crop, location)) * (1 - exp(-irrigParam(crop, location)*irrigI(crop, location))) +
+ (yieldBoth(crop, location) + yieldNone(crop, location) - yieldFertOnly(crop, location) - yieldIrrigOnly(crop, location)) *
+ (1 - exp(-fertParam(crop, location)*fertI(crop, location))) * (1 - exp(-irrigParam(crop, location)*irrigI(crop, location)))
+ ) * (1 - exp(-otherIntensity(crop, location)*otherIParam));
+
+ NET_SUPPLY_EQ(crop) .. net_supply(crop) =E= sum(location, cropArea(crop, location) * yield(crop, location)) * (1 - seedAndWasteRate(crop)) - ruminantFeed(crop) - monogastricFeed(crop) + importAmount(crop) - exportAmount(crop);
+
+ CROP_DEMAND_CONSTRAINT(crop) .. net_supply(crop) =G= demand(crop);
+
+ TOTAL_CEREAL_DEMAND_CONSTRAINT .. sum(cereal_crop, net_supply(cereal_crop)) =G= demand('cereals');
+
+ TOTAL_OIL_PULSE_DEMAND_CONSTRAINT .. sum(oilpulse_crop, net_supply(oilpulse_crop)) =G= demand('oilcropspulses');
+
+ RUMINANT_DEMAND_CONSTRAINT .. meatEfficency * sum(feed_crop, ruminantFeed(feed_crop) * cropDM(feed_crop)) * (1 - seedAndWasteRate('ruminants')) =G= (demand('ruminants') - importAmount('ruminants') + exportAmount('ruminants'));
+
+ MONOGASTRICS_DEMAND_CONSTRAINT .. meatEfficency * sum(feed_crop_less_pasture, monogastricFeed(feed_crop_less_pasture) * cropDM(feed_crop_less_pasture)) * (1 - seedAndWasteRate('monogastrics')) =G= (demand('monogastrics') - importAmount('monogastrics') + exportAmount('monogastrics'));
+
+ TOTAL_NON_PASTURE_FEED_DM_CALC .. totalFeedDM =E= sum(feed_crop_less_pasture, (ruminantFeed(feed_crop_less_pasture) + monogastricFeed(feed_crop_less_pasture)) * cropDM(feed_crop_less_pasture));
+ FEED_MIX_CONSTRAINT(feed_crop_less_pasture) .. (ruminantFeed(feed_crop_less_pasture) + monogastricFeed(feed_crop_less_pasture)) * cropDM(feed_crop_less_pasture) =L= totalFeedDM * 0.7;
+
+ MAX_FERT_INTENSITY_CONSTRAINT(crop, location) .. fertI(crop, location) =L= 1;
+ MAX_IRRIG_INTENSITY_CONSTRAINT(crop, location) .. irrigI(crop, location) =L= 1;
+ MAX_OTHER_INTENSITY_CONSTRAINT(crop, location) .. otherIntensity(crop, location) =L= 1;
+
+ SETASIDE_AREA_CALC(location) .. cropArea('setaside', location) =E= sum(crop_less_pasture, cropArea(crop_less_pasture, location)) * setAsideRate;
+
+ MAX_NET_IMPORT_CONSTRAINT(import_crop) .. importAmount(import_crop) - exportAmount(import_crop) =L= maxNetImport(import_crop);
+ MIN_NET_IMPORT_CONSTRAINT(import_crop) .. importAmount(import_crop) - exportAmount(import_crop) =G= minNetImport(import_crop);
+ PASTURE_IMPORT_CONSTRAINT .. importAmount('pasture') =E= 0;
+ PASTURE_EXPORT_CONSTRAINT .. exportAmount('pasture') =E= 0;
+
+ IRRIGATION_CONSTRAINT(location) .. irrigConstraint(location) * suitableLandArea(location) * (1.0 - unhandledCropRate) =G= sum(crop, irrigMaxRate(crop, location) * irrigI(crop, location) * cropArea(crop, location));
+
+ CROPLAND_LAND_COVER_CALC(location) .. landCoverArea('cropland', location) =E= sum(crop_less_pasture, cropArea(crop_less_pasture, location)) / (1.0 - unhandledCropRate);
+
+ PASTURE_LAND_COVER_CALC(location) .. landCoverArea('pasture', location) =E= cropArea('pasture', 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));
+
+* landCoverChange(A, A, location) defined as unchanged land cover
+ LAND_COVER_CHANGE_CONSTRAINT(land_cover, location) .. sum(land_cover_after, landCoverChange(land_cover, land_cover_after, location)) =E= previousLandCoverArea(land_cover, location);
+
+ NATURAL_CONVERSION_CONSTRAINT(exc_natural, location) .. landCoverChange(exc_natural, 'natural', location) =E= 0;
+
+ NATURAL_CONTRACTION_CONSTRAINT(location) .. landCoverArea('natural', location) =G= previousLandCoverArea('natural', location) * (1 - naturalDecreaseMaxRate);
+
+ CONVERSION_COST(location) .. totalConversionCost(location) =E= sum((land_cover_before, land_cover_after), landCoverChange(land_cover_before, land_cover_after, location) * conversionCost(land_cover_before, land_cover_after)) +
+ sum(exc_natural, landCoverChange('natural', exc_natural, location) * infrastructureCost(location));
+
+* Wood from land cover change
+* WOOD_HARVEST_LUC_CALC(location) .. woodHarvestLuc(location) =E= sum((land_cover_before, land_cover_after), landCoverChange(land_cover_before, land_cover_after, location) * woodYieldLuc(land_cover_before, land_cover_after, location));
+
+ WOOD_HARVEST_LUC_CALC(location) .. woodHarvestLuc(location) =E= woodYieldLuc('natural', 'naturalNew', location) * landCoverChange('natural', 'naturalNew', location);
+
+* Wood from timberForest rotation = harvest from new forest + existing forest
+ WOOD_HARVEST_ROTA_CALC(location) .. woodHarvestRota(location) =E= [ sum(non_timber_forest, landCoverChange(non_timber_forest, 'timberForest', location) * woodYieldRota(non_timber_forest, 'timberForest', location)) +
+ landCoverChange('timberForest', 'timberForest', location) * timberForestYield(location)
+ ] / forestRotationPeriod;
+
+ TOTAL_WOOD_HARVEST_CALC(location) .. totalWoodHarvest(location) =E= woodHarvestLuc(location) + woodHarvestRota(location);
+
+ NATURAL_WOOD_HARVEST_CONSTRAINT(location) .. woodHarvestLuc(location) =L= previousLandCoverArea('natural', location) * woodYieldLuc('natural', 'naturalNew', location) * naturalWoodHarvestMaxRate;
+
+ WOOD_MAX_EXPORT_CONSTRAINT .. woodExported =L= sum(location, totalWoodHarvest(location));
+
+ WOOD_MIN_TRADE_CONSTRAINT .. 0 - woodExported =L= timberMaxNetImport;
+
+ WOOD_MAX_TRADE_CONSTRAINT .. 0 - woodExported =G= timberMinNetImport;
+
+* Must sell timber producted in timberForest
+ WOOD_ROTA_EXPORT_CONSTRAINT .. woodExported =G= sum(location, woodHarvestRota(location));
+
+* WOOD_ROTA_TO_LUC_CONSTRAINT .. sum(location, woodHarvestLuc(location)) =L= sum(location, woodHarvestRota(location)) * 2;
+
+ LOGGING_COST_CALC(location) .. loggingCost(location) =E= 0.04 + 0.06 * landCoverChange('natural', 'naturalNew', location) / [previousLandCoverArea('natural', location) + delta] *
+ landCoverChange('natural', 'naturalNew', location) * woodYieldLuc('natural', 'naturalNew', location);
+
+ CARBON_FLUX_CALC(location) .. carbonFlux(location) =E= sum((land_cover_before, land_cover_after), landCoverChange(land_cover_before, land_cover_after, location) * carbonFluxRate(land_cover_before, land_cover_after, location));
+
+
+ COST_EQ .. total_cost =E=
+ (
+ ( SUM((crop, location), cropArea(crop, location) * unitCost(crop, location) * (1-subsidyRate(crop))) ) * domesticPriceMarkup +
+
+ SUM(import_crop, importAmount(import_crop) * importPrices(import_crop) - exportAmount(import_crop) * exportPrices(import_crop)) +
+
+ SUM(location, totalConversionCost(location)) +
+
+ [SUM((managed_forest, location), landCoverArea(managed_forest, location) * forestRotationCost) - woodExported * woodPrice] + SUM(location, loggingCost(location)) +
+
+ SUM(location, carbonFlux(location)) * carbonPrice
+ );
+
+ MODEL LAND_USE /ALL/ ;
+ fertI.L(crop, location) = previousFertIntensity(crop, location);
+ irrigI.L(crop, location) = previousIrrigIntensity(crop, location);
+ otherIntensity.L(crop, location) = previousOtherIntensity(crop, location);
+ cropArea.L(crop, location) = previousCropArea(crop, location);
+ landCoverArea.L(land_cover, location) = previousLandCoverArea(land_cover, location);
+ landCoverChange.L(land_cover, land_cover, location) = previousLandCoverArea(land_cover, location);
+ ruminantFeed.L(feed_crop) = previousRuminantFeed(feed_crop);
+ monogastricFeed.L(feed_crop) = previousMonogastricFeed(feed_crop);
+ importAmount.L(all_types) = previousImportAmount(all_types);
+ exportAmount.L(all_types) = previousExportAmount(all_types);
+ woodHarvestRota.L(location) = previousLandCoverArea('timberForest', location) * timberForestYield(location) / forestRotationPeriod;
+ totalConversionCost.L(location) = previousLandCoverArea('timberForest', location) * conversionCost('timberForest', 'timberForest');
+
+ LAND_USE.OptFile = 1;
+
+* display landCoverChange.L
+
+ SOLVE LAND_USE USING NLP MINIMIZING total_cost;
+
+ display previousCropArea, irrigMaxRate, otherIntensity.L, fertI.L, irrigI.L, cropArea.L;
+ display net_supply.l, demand, ruminantFeed.l, monogastricFeed.l, importAmount.l, exportAmount.l;
+
+* Calculate summary information used in Java process
+ parameter totalProd(all_types);
+ parameter totalProdCost(all_types);
+ parameter totalArea(crop);
+ parameter totalCropland(location);
+ parameter netImportAmount(all_types);
+ parameter netImportCost(all_types);
+ parameter feedCostRate(feed_crop);
+ parameter productionShock(all_types);
+ scalar netCarbonFlux;
+
+* Production quantities based on smaller area (before unhandledCropArea adjustment applied)
+ totalProd(crop) = sum(location, cropArea.l(crop, location) * yield.l(crop, location));
+ productionShock(crop) = sum(location, cropArea.l(crop, location) * yield.l(crop, location) * yieldShock(crop, location));
+ totalProd('ruminants') = meatEfficency*(sum(feed_crop, ruminantFeed.l(feed_crop) * cropDM(feed_crop)));
+ totalProd('monogastrics') = meatEfficency*(sum(feed_crop, monogastricFeed.l(feed_crop) * cropDM(feed_crop)));
+
+* Cost based on adjusted area
+ cropArea.l(crop_less_pasture, location) = cropArea.l(crop_less_pasture, location) / (1.0 - unhandledCropRate);
+ totalProdCost(crop) = sum(location, unitCost.l(crop, location) * cropArea.l(crop, location));
+ totalCropland(location) = sum(crop_less_pasture, cropArea.l(crop_less_pasture, location));
+ totalArea(crop) = sum(location, cropArea.l(crop, location));
+
+ feedCostRate(feed_crop)$[totalProd(feed_crop) <> 0] = totalProdCost(feed_crop) / totalProd(feed_crop);
+ totalProdCost('ruminants') = sum(feed_crop, ruminantFeed.l(feed_crop) * feedCostRate(feed_crop));
+ totalProdCost('monogastrics') = sum(feed_crop, monogastricFeed.l(feed_crop) * feedCostRate(feed_crop));
+
+ netImportCost(import_crop) = importAmount.l(import_crop) * importPrices(import_crop) - exportAmount.l(import_crop) * exportPrices(import_crop);
+ netImportAmount(import_crop) = importAmount.l(import_crop) - exportAmount.l(import_crop);
+
+ netCarbonFlux = SUM(location, carbonFlux.L(location));
+
+ Scalar totalCostsLessLU;
+
+ totalCostsLessLU = sum(crop,totalProdCost(crop))+ sum(import_crop,netImportCost(import_crop));
+
+
+ Scalar ms 'model status', ss 'solve status';
+ ms=LAND_USE.modelstat;
+ ss=LAND_USE.solvestat;
diff --git a/data/timber_demand.csv b/data/timber_demand.csv
index 467ef612314c2da2af6642371479f64a09f8e67e..d2d6519f06f221b9d5e591f75b2f7bc4325a7db7 100644
--- a/data/timber_demand.csv
+++ b/data/timber_demand.csv
@@ -1,22 +1,387 @@
-year,value
-2000,945
-2005,954
-2010,964
-2015,973
-2020,983
-2025,993
-2030,1003
-2035,1013
-2040,1023
-2045,1033
-2050,1043
-2055,1054
-2060,1064
-2065,1075
-2070,1086
-2075,1097
-2080,1108
-2085,1119
-2090,1130
-2095,1141
-2100,1153
+country,type,demand m3
+New Zealand,roundwood,13747205
+Australia,roundwood,24114981
+New Caledonia,roundwood,15289
+China,roundwood,197175942
+Tuvalu,roundwood,64
+Palau,roundwood,632
+Marshall Islands,roundwood,0
+Kiribati,roundwood,6
+Micronesia (Federated States of),roundwood,493
+Tonga,roundwood,16586
+Samoa,roundwood,6317
+Vanuatu,roundwood,37898
+Brunei Darussalam,roundwood,107187
+Solomon Islands,roundwood,30000
+Macau,roundwood,0
+Fiji,roundwood,796286
+Lao People's Democratic Republic,roundwood,825266
+Papua New Guinea,roundwood,1691789
+French Polynesia,roundwood,3466
+Democratic People's Republic of Korea,roundwood,1490232
+Indonesia,roundwood,54131388
+Japan,roundwood,21884141
+Mongolia,roundwood,48833
+Myanmar,roundwood,2306527
+Singapore,roundwood,13648
+Malaysia,roundwood,14846894
+Philippines,roundwood,3599731
+Republic of Korea,roundwood,7682740
+Thailand,roundwood,14859173
+Cambodia,roundwood,282362
+Viet Nam,roundwood,24982678
+Albania,roundwood,81393
+Armenia,roundwood,987
+Slovakia,roundwood,7237478
+Serbia,roundwood,1485310
+Bulgaria,roundwood,2563824
+Hungary,roundwood,2134211
+Czechia,roundwood,14007000
+Romania,roundwood,10790478
+Montenegro,roundwood,194341
+Kosovo,roundwood,0
+Slovenia,roundwood,1469301
+The former Yugoslav Republic of Macedonia,roundwood,134814
+Republic of Moldova,roundwood,79452
+Bosnia and Herzegovina,roundwood,2355090
+Croatia,roundwood,2849000
+Luxembourg,roundwood,1111821
+Ireland,roundwood,2204027
+Belgium,roundwood,7451697
+Netherlands,roundwood,519793
+France,roundwood,24659540
+Switzerland,roundwood,2931741
+Austria,roundwood,20367313
+Germany,roundwood,51066044
+Italy,roundwood,5798688
+Georgia,roundwood,120800
+Turkmenistan,roundwood,0
+Kyrgyzstan,roundwood,12600
+Tajikistan,roundwood,0
+Azerbaijan,roundwood,13990
+Belarus,roundwood,5890206
+Kazakhstan,roundwood,158136
+Uzbekistan,roundwood,278257
+Poland,roundwood,32047012
+Estonia,roundwood,2891583
+Latvia,roundwood,6500644
+Lithuania,roundwood,4105388
+Russian Federation,roundwood,140613360
+Greenland,roundwood,72
+Iceland,roundwood,150
+Norway,roundwood,8745709
+Finland,roundwood,51192349
+Denmark,roundwood,1847218
+Sweden,roundwood,71359465
+Portugal,roundwood,8902821
+Spain,roundwood,11476645
+Cyprus,roundwood,5687
+Greece,roundwood,468174
+Turkey,roundwood,16861600
+Ukraine,roundwood,4622300
+United Kingdom,roundwood,8215840
+Argentina,roundwood,10589282
+Bolivia (Plurinational State of),roundwood,937091
+Chile,roundwood,34525165
+Brazil,roundwood,128399924
+Saint Kitts and Nevis,roundwood,351
+Dominica,roundwood,831
+Antigua and Barbuda,roundwood,377
+Aruba,roundwood,2667
+Grenada,roundwood,0
+Saint Vincent and the Grenadines,roundwood,5507
+Saint Lucia,roundwood,925
+Barbados,roundwood,7158
+Bahamas,roundwood,24096
+Trinidad and Tobago,roundwood,65170
+Jamaica,roundwood,166939
+Puerto Rico,roundwood,0
+Dominican Republic,roundwood,22287
+Haiti,roundwood,263186
+Cuba,roundwood,605642
+Bermuda,roundwood,0
+Belize,roundwood,46542
+Panama,roundwood,132349
+Costa Rica,roundwood,1212485
+Nicaragua,roundwood,104707
+El Salvador,roundwood,658854
+Honduras,roundwood,618142
+Guatemala,roundwood,682453
+Colombia,roundwood,3532425
+Mexico,roundwood,5054930
+Uruguay,roundwood,8936000
+Paraguay,roundwood,4027946
+Ecuador,roundwood,1980558
+Peru,roundwood,1352355
+Suriname,roundwood,197830
+Guyana,roundwood,400426
+Venezuela (Bolivarian Republic of),roundwood,1850401
+Algeria,roundwood,154686
+Egypt,roundwood,399388
+Iran (Islamic Republic of),roundwood,713945
+Iraq,roundwood,60353
+Bahrain,roundwood,3995
+Qatar,roundwood,3543
+Oman,roundwood,12331
+Kuwait,roundwood,107
+Lebanon,roundwood,25986
+Jordan,roundwood,5389
+Israel,roundwood,29653
+United Arab Emirates,roundwood,158934
+Yemen,roundwood,1271
+Saudi Arabia,roundwood,23430
+Morocco,roundwood,610270
+Malta,roundwood,296
+Djibouti,roundwood,1416
+Libya,roundwood,132306
+Tunisia,roundwood,232586
+Canada,roundwood,139528000
+United States of America,roundwood,325368790
+Afghanistan,roundwood,1759661
+Bangladesh,roundwood,436733
+Maldives,roundwood,44
+Sri Lanka,roundwood,698594
+India,roundwood,54055057
+Bhutan,roundwood,189346
+Nepal,roundwood,1299992
+Pakistan,roundwood,4201869
+Sao Tome and Principe,roundwood,27120
+Comoros,roundwood,24650
+Equatorial Guinea,roundwood,83038
+Congo,roundwood,1429442
+Central African Republic,roundwood,484390
+Burundi,roundwood,880507
+Rwanda,roundwood,1213085
+Chad,roundwood,761508
+Cameroon,roundwood,1646546
+Democratic Republic of the Congo,roundwood,4369117
+Eritrea,roundwood,5758
+South Sudan,roundwood,540882
+Ethiopia,roundwood,2936000
+Kenya,roundwood,1147092
+Nigeria,roundwood,9308104
+South Africa,roundwood,16818652
+Swaziland,roundwood,935117
+Mauritius,roundwood,32051
+Lesotho,roundwood,2995
+Botswana,roundwood,113360
+Namibia,roundwood,4266
+Zambia,roundwood,1451780
+Zimbabwe,roundwood,518434
+Malawi,roundwood,1474352
+Madagascar,roundwood,284002
+Angola,roundwood,1124316
+Mozambique,roundwood,1180044
+Sudan,roundwood,1622648
+Uganda,roundwood,4093645
+United Republic of Tanzania,roundwood,2305005
+Cabo Verde,roundwood,1274
+Gabon,roundwood,895000
+Guinea-Bissau,roundwood,116947
+Gambia,roundwood,113350
+Mauritania,roundwood,50058
+Liberia,roundwood,472271
+Sierra Leone,roundwood,108715
+Togo,roundwood,48579
+Benin,roundwood,380600
+Guinea,roundwood,501817
+Senegal,roundwood,791761
+Mali,roundwood,412549
+Burkina Faso,roundwood,1170463
+Niger,roundwood,701543
+Cote d'Ivoire,roundwood,2210854
+Ghana,roundwood,1808723
+Somalia,roundwood,111187
+New Zealand,fuelwood,-96
+Australia,fuelwood,4852740
+New Caledonia,fuelwood,11819
+China,fuelwood,188825892
+Tuvalu,fuelwood,0
+Palau,fuelwood,0
+Marshall Islands,fuelwood,0
+Kiribati,fuelwood,3058
+Micronesia (Federated States of),fuelwood,2659
+Tonga,fuelwood,2137
+Samoa,fuelwood,70014
+Vanuatu,fuelwood,90404
+Brunei Darussalam,fuelwood,11747
+Solomon Islands,fuelwood,128451
+Macau,fuelwood,0
+Fiji,fuelwood,36989
+Lao People's Democratic Republic,fuelwood,5947762
+Papua New Guinea,fuelwood,5533000
+French Polynesia,fuelwood,4513
+Democratic People's Republic of Korea,fuelwood,5986742
+Indonesia,fuelwood,59796369
+Japan,fuelwood,87567
+Mongolia,fuelwood,763403
+Myanmar,fuelwood,38285981
+Singapore,fuelwood,32613
+Malaysia,fuelwood,2805619
+Philippines,fuelwood,12361935
+Republic of Korea,fuelwood,518047
+Thailand,fuelwood,19307543
+Cambodia,fuelwood,8442343
+Viet Nam,fuelwood,21467259
+Albania,fuelwood,293700
+Armenia,fuelwood,1750000
+Slovakia,fuelwood,448124
+Serbia,fuelwood,6213400
+Bulgaria,fuelwood,2464625
+Hungary,fuelwood,2996397
+Czechia,fuelwood,1939000
+Romania,fuelwood,2502956
+Montenegro,fuelwood,704087
+Kosovo,fuelwood,0
+Slovenia,fuelwood,938993
+The former Yugoslav Republic of Macedonia,fuelwood,572935
+Republic of Moldova,fuelwood,311164
+Bosnia and Herzegovina,fuelwood,929830
+Croatia,fuelwood,816000
+Luxembourg,fuelwood,78659
+Ireland,fuelwood,185287
+Belgium,fuelwood,758013
+Netherlands,fuelwood,288000
+France,fuelwood,25403035
+Switzerland,fuelwood,1485022
+Austria,fuelwood,5084769
+Germany,fuelwood,27578693
+Italy,fuelwood,12380000
+Georgia,fuelwood,733000
+Turkmenistan,fuelwood,10000
+Kyrgyzstan,fuelwood,36449
+Tajikistan,fuelwood,90000
+Azerbaijan,fuelwood,3200
+Belarus,fuelwood,2287481
+Kazakhstan,fuelwood,272070
+Uzbekistan,fuelwood,22000
+Poland,fuelwood,4010696
+Estonia,fuelwood,2028855
+Latvia,fuelwood,984536
+Lithuania,fuelwood,1881659
+Russian Federation,fuelwood,13710707
+Greenland,fuelwood,9
+Iceland,fuelwood,34
+Norway,fuelwood,2292316
+Finland,fuelwood,6838296
+Denmark,fuelwood,1498661
+Sweden,fuelwood,6375717
+Portugal,fuelwood,598838
+Spain,fuelwood,5062690
+Cyprus,fuelwood,3865
+Greece,fuelwood,1007658
+Turkey,fuelwood,5143988
+Ukraine,fuelwood,7871741
+United Kingdom,fuelwood,1225345
+Argentina,fuelwood,4374361
+Bolivia (Plurinational State of),fuelwood,2349841
+Chile,fuelwood,12655000
+Brazil,fuelwood,107032114
+Saint Kitts and Nevis,fuelwood,0
+Dominica,fuelwood,7493
+Antigua and Barbuda,fuelwood,0
+Aruba,fuelwood,1749
+Grenada,fuelwood,0
+Saint Vincent and the Grenadines,fuelwood,7444
+Saint Lucia,fuelwood,9850
+Barbados,fuelwood,4982
+Bahamas,fuelwood,32907
+Trinidad and Tobago,fuelwood,32875
+Jamaica,fuelwood,545053
+Puerto Rico,fuelwood,0
+Dominican Republic,fuelwood,912647
+Haiti,fuelwood,2041142
+Cuba,fuelwood,1141000
+Bermuda,fuelwood,0
+Belize,fuelwood,126079
+Panama,fuelwood,1128086
+Costa Rica,fuelwood,3377015
+Nicaragua,fuelwood,6095779
+El Salvador,fuelwood,4225393
+Honduras,fuelwood,8574891
+Guatemala,fuelwood,18058663
+Colombia,fuelwood,8826028
+Mexico,fuelwood,38820610
+Uruguay,fuelwood,2430000
+Paraguay,fuelwood,6575817
+Ecuador,fuelwood,4940180
+Peru,fuelwood,7338000
+Suriname,fuelwood,96704
+Guyana,fuelwood,848053
+Venezuela (Bolivarian Republic of),fuelwood,4054562
+Algeria,fuelwood,8176482
+Egypt,fuelwood,17511977
+Iran (Islamic Republic of),fuelwood,53972
+Iraq,fuelwood,118000
+Bahrain,fuelwood,6844
+Qatar,fuelwood,5348
+Oman,fuelwood,38943
+Kuwait,fuelwood,18504
+Lebanon,fuelwood,20460
+Jordan,fuelwood,302609
+Israel,fuelwood,2130
+United Arab Emirates,fuelwood,22624
+Yemen,fuelwood,440876
+Saudi Arabia,fuelwood,249039
+Morocco,fuelwood,6785858
+Malta,fuelwood,25
+Djibouti,fuelwood,355976
+Libya,fuelwood,951892
+Tunisia,fuelwood,2184631
+Canada,fuelwood,3255879
+United States of America,fuelwood,40262948
+Afghanistan,fuelwood,1641303
+Bangladesh,fuelwood,27286998
+Maldives,fuelwood,15577
+Sri Lanka,fuelwood,5211428
+India,fuelwood,309313649
+Bhutan,fuelwood,4845000
+Nepal,fuelwood,12525698
+Pakistan,fuelwood,29659986
+Sao Tome and Principe,fuelwood,107428
+Comoros,fuelwood,265913
+Equatorial Guinea,fuelwood,447000
+Congo,fuelwood,1335790
+Central African Republic,fuelwood,2000000
+Burundi,fuelwood,9259281
+Rwanda,fuelwood,4999980
+Chad,fuelwood,7070034
+Cameroon,fuelwood,9905899
+Democratic Republic of the Congo,fuelwood,76602014
+Eritrea,fuelwood,920987
+South Sudan,fuelwood,18776039
+Ethiopia,fuelwood,101274510
+Kenya,fuelwood,26402590
+Nigeria,fuelwood,63215134
+South Africa,fuelwood,12268000
+Swaziland,fuelwood,1062471
+Mauritius,fuelwood,8997
+Lesotho,fuelwood,2092623
+Botswana,fuelwood,683425
+Namibia,fuelwood,1307674
+Zambia,fuelwood,21822922
+Zimbabwe,fuelwood,8708935
+Malawi,fuelwood,5404338
+Madagascar,fuelwood,13100000
+Angola,fuelwood,4009343
+Mozambique,fuelwood,16724246
+Sudan,fuelwood,18776039
+Uganda,fuelwood,39636345
+United Republic of Tanzania,fuelwood,22833027
+Cabo Verde,fuelwood,186055
+Gabon,fuelwood,1069998
+Guinea-Bissau,fuelwood,2599954
+Gambia,fuelwood,694278
+Mauritania,fuelwood,1836252
+Liberia,fuelwood,7008341
+Sierra Leone,fuelwood,5581587
+Togo,fuelwood,4423968
+Benin,fuelwood,6274504
+Guinea,fuelwood,11958951
+Senegal,fuelwood,5427207
+Mali,fuelwood,5325843
+Burkina Faso,fuelwood,12784965
+Niger,fuelwood,9875857
+Cote d'Ivoire,fuelwood,8945530
+Ghana,fuelwood,37782939
+Somalia,fuelwood,13500863
diff --git a/data/timber_demand_world.csv b/data/timber_demand_world.csv
index 467ef612314c2da2af6642371479f64a09f8e67e..7735b564968bc750badc179d931680a295cdc61e 100644
--- a/data/timber_demand_world.csv
+++ b/data/timber_demand_world.csv
@@ -1,22 +1,387 @@
-year,value
-2000,945
-2005,954
-2010,964
-2015,973
-2020,983
-2025,993
-2030,1003
-2035,1013
-2040,1023
-2045,1033
-2050,1043
-2055,1054
-2060,1064
-2065,1075
-2070,1086
-2075,1097
-2080,1108
-2085,1119
-2090,1130
-2095,1141
-2100,1153
+country,type,demand m3
+New Zealand,roundwood,24490000
+Australia,roundwood,25561000
+New Caledonia,roundwood,14700
+China,roundwood,161809500
+Tuvalu,roundwood,0
+Palau,roundwood,0
+Marshall Islands,roundwood,0
+Kiribati,roundwood,0
+Micronesia (Federated States of),roundwood,0
+Tonga,roundwood,2000
+Samoa,roundwood,5700
+Vanuatu,roundwood,38000
+Brunei Darussalam,roundwood,107300
+Solomon Islands,roundwood,1624000
+Macau,roundwood,0
+Fiji,roundwood,800000
+Lao People's Democratic Republic,roundwood,1032000
+Papua New Guinea,roundwood,4476000
+French Polynesia,roundwood,1000
+Democratic People's Republic of Korea,roundwood,1500000
+Indonesia,roundwood,54105500
+Japan,roundwood,17193000
+Mongolia,roundwood,49000
+Myanmar,roundwood,4171000
+Singapore,roundwood,0
+Malaysia,roundwood,19303000
+Philippines,roundwood,3593000
+Republic of Korea,roundwood,3466000
+Thailand,roundwood,14600000
+Cambodia,roundwood,288000
+Viet Nam,roundwood,24200000
+Albania,roundwood,80000
+Armenia,roundwood,1000
+Slovakia,roundwood,9089175
+Serbia,roundwood,1413000
+Bulgaria,roundwood,3011000
+Hungary,roundwood,2746292
+Czechia,roundwood,14771000
+Romania,roundwood,10548052
+Montenegro,roundwood,208001
+Kosovo,roundwood,0
+Slovenia,roundwood,1841404
+The former Yugoslav Republic of Macedonia,roundwood,101000
+Republic of Moldova,roundwood,43000
+Bosnia and Herzegovina,roundwood,2461000
+Croatia,roundwood,3421000
+Luxembourg,roundwood,420593
+Ireland,roundwood,2436975
+Belgium,roundwood,4113900
+Netherlands,roundwood,790593
+France,roundwood,29634241
+Switzerland,roundwood,3439389
+Austria,roundwood,13281444
+Germany,roundwood,47136000
+Italy,roundwood,2647228
+Georgia,roundwood,105006
+Turkmenistan,roundwood,0
+Kyrgyzstan,roundwood,9300
+Tajikistan,roundwood,0
+Azerbaijan,roundwood,3300
+Belarus,roundwood,8072600
+Kazakhstan,roundwood,73000
+Uzbekistan,roundwood,8000
+Poland,roundwood,31343001
+Estonia,roundwood,4841309
+Latvia,roundwood,10221818
+Lithuania,roundwood,5153862
+Russian Federation,roundwood,161595000
+Greenland,roundwood,0
+Iceland,roundwood,0
+Norway,roundwood,8322430
+Finland,roundwood,45419600
+Denmark,roundwood,1695300
+Sweden,roundwood,66300000
+Portugal,roundwood,9048360
+Spain,roundwood,10969399
+Cyprus,roundwood,5330
+Greece,roundwood,336481
+Turkey,roundwood,15695000
+Ukraine,roundwood,7536000
+United Kingdom,roundwood,8337163
+Argentina,roundwood,10635024
+Bolivia (Plurinational State of),roundwood,943000
+Chile,roundwood,34560000
+Brazil,roundwood,128400000
+Saint Kitts and Nevis,roundwood,0
+Dominica,roundwood,0
+Antigua and Barbuda,roundwood,0
+Aruba,roundwood,0
+Grenada,roundwood,0
+Saint Vincent and the Grenadines,roundwood,0
+Saint Lucia,roundwood,0
+Barbados,roundwood,6000
+Bahamas,roundwood,17000
+Trinidad and Tobago,roundwood,65000
+Jamaica,roundwood,154805
+Puerto Rico,roundwood,0
+Dominican Republic,roundwood,9700
+Haiti,roundwood,239000
+Cuba,roundwood,589500
+Bermuda,roundwood,0
+Belize,roundwood,41000
+Panama,roundwood,176000
+Costa Rica,roundwood,1326000
+Nicaragua,roundwood,108000
+El Salvador,roundwood,682000
+Honduras,roundwood,620000
+Guatemala,roundwood,692000
+Colombia,roundwood,3550000
+Mexico,roundwood,5041000
+Uruguay,roundwood,9402000
+Paraguay,roundwood,4044000
+Ecuador,roundwood,2091000
+Peru,roundwood,1352000
+Suriname,roundwood,246300
+Guyana,roundwood,516000
+Venezuela (Bolivarian Republic of),roundwood,1850000
+Algeria,roundwood,144600
+Egypt,roundwood,268000
+Iran (Islamic Republic of),roundwood,697000
+Iraq,roundwood,59000
+Bahrain,roundwood,0
+Qatar,roundwood,0
+Oman,roundwood,0
+Kuwait,roundwood,0
+Lebanon,roundwood,18000
+Jordan,roundwood,4000
+Israel,roundwood,24957
+United Arab Emirates,roundwood,0
+Yemen,roundwood,0
+Saudi Arabia,roundwood,0
+Morocco,roundwood,372100
+Malta,roundwood,0
+Djibouti,roundwood,0
+Libya,roundwood,116000
+Tunisia,roundwood,218000
+Canada,roundwood,138802000
+United States of America,roundwood,336135000
+Afghanistan,roundwood,1760000
+Bangladesh,roundwood,420000
+Maldives,roundwood,0
+Sri Lanka,roundwood,698000
+India,roundwood,48759000
+Bhutan,roundwood,190000
+Nepal,roundwood,1300000
+Pakistan,roundwood,4060000
+Sao Tome and Principe,roundwood,27200
+Comoros,roundwood,24650
+Equatorial Guinea,roundwood,309849
+Congo,roundwood,2035524
+Central African Republic,roundwood,632283
+Burundi,roundwood,883000
+Rwanda,roundwood,1211927
+Chad,roundwood,761000
+Cameroon,roundwood,2348000
+Democratic Republic of the Congo,roundwood,4531539
+Eritrea,roundwood,5654
+South Sudan,roundwood,0
+Ethiopia,roundwood,2935000
+Kenya,roundwood,1133000
+Nigeria,roundwood,9418000
+South Africa,roundwood,16988569
+Swaziland,roundwood,934000
+Mauritius,roundwood,5000
+Lesotho,roundwood,0
+Botswana,roundwood,105000
+Namibia,roundwood,0
+Zambia,roundwood,1455000
+Zimbabwe,roundwood,518400
+Malawi,roundwood,1480000
+Madagascar,roundwood,280911
+Angola,roundwood,1127000
+Mozambique,roundwood,1416000
+Sudan,roundwood,2173000
+Uganda,roundwood,4093000
+United Republic of Tanzania,roundwood,2314000
+Cabo Verde,roundwood,0
+Gabon,roundwood,1800000
+Guinea-Bissau,roundwood,131910
+Gambia,roundwood,131700
+Mauritania,roundwood,32000
+Liberia,roundwood,480000
+Sierra Leone,roundwood,123600
+Togo,roundwood,166000
+Benin,roundwood,467000
+Guinea,roundwood,651000
+Senegal,roundwood,788000
+Mali,roundwood,412900
+Burkina Faso,roundwood,1171000
+Niger,roundwood,701000
+Cote d'Ivoire,roundwood,2356000
+Ghana,roundwood,1950000
+Somalia,roundwood,110000
+New Zealand,fuelwood,0
+Australia,fuelwood,4853000
+New Caledonia,fuelwood,11818
+China,fuelwood,188823052
+Tuvalu,fuelwood,0
+Palau,fuelwood,0
+Marshall Islands,fuelwood,0
+Kiribati,fuelwood,2937
+Micronesia (Federated States of),fuelwood,2444
+Tonga,fuelwood,2134
+Samoa,fuelwood,70000
+Vanuatu,fuelwood,91000
+Brunei Darussalam,fuelwood,11679
+Solomon Islands,fuelwood,128451
+Macau,fuelwood,0
+Fiji,fuelwood,37000
+Lao People's Democratic Republic,fuelwood,5947951
+Papua New Guinea,fuelwood,5533000
+French Polynesia,fuelwood,4417
+Democratic People's Republic of Korea,fuelwood,5986734
+Indonesia,fuelwood,59743242
+Japan,fuelwood,87729
+Mongolia,fuelwood,763361
+Myanmar,fuelwood,38286000
+Singapore,fuelwood,32799
+Malaysia,fuelwood,2810489
+Philippines,fuelwood,12361935
+Republic of Korea,fuelwood,518000
+Thailand,fuelwood,19301440
+Cambodia,fuelwood,8442343
+Viet Nam,fuelwood,21500000
+Albania,fuelwood,350000
+Armenia,fuelwood,1750000
+Slovakia,fuelwood,509893
+Serbia,fuelwood,6223000
+Bulgaria,fuelwood,2657000
+Hungary,fuelwood,2993983
+Czechia,fuelwood,1965000
+Romania,fuelwood,2563588
+Montenegro,fuelwood,707000
+Kosovo,fuelwood,0
+Slovenia,fuelwood,1104045
+The former Yugoslav Republic of Macedonia,fuelwood,530000
+Republic of Moldova,fuelwood,308800
+Bosnia and Herzegovina,fuelwood,1418000
+Croatia,fuelwood,1056000
+Luxembourg,fuelwood,81576
+Ireland,fuelwood,181021
+Belgium,fuelwood,713525
+Netherlands,fuelwood,290000
+France,fuelwood,26173567
+Switzerland,fuelwood,1498647
+Austria,fuelwood,4549512
+Germany,fuelwood,27296000
+Italy,fuelwood,11429000
+Georgia,fuelwood,733000
+Turkmenistan,fuelwood,10000
+Kyrgyzstan,fuelwood,36600
+Tajikistan,fuelwood,90000
+Azerbaijan,fuelwood,3200
+Belarus,fuelwood,2291600
+Kazakhstan,fuelwood,272000
+Uzbekistan,fuelwood,22000
+Poland,fuelwood,4124416
+Estonia,fuelwood,2195874
+Latvia,fuelwood,2312000
+Lithuania,fuelwood,1942998
+Russian Federation,fuelwood,13904000
+Greenland,fuelwood,0
+Iceland,fuelwood,0
+Norway,fuelwood,2120652
+Finland,fuelwood,6705000
+Denmark,fuelwood,1362200
+Sweden,fuelwood,5900000
+Portugal,fuelwood,600000
+Spain,fuelwood,5120000
+Cyprus,fuelwood,3628
+Greece,fuelwood,711481
+Turkey,fuelwood,4902000
+Ukraine,fuelwood,8609600
+United Kingdom,fuelwood,1381100
+Argentina,fuelwood,4374522
+Bolivia (Plurinational State of),fuelwood,2349841
+Chile,fuelwood,12655000
+Brazil,fuelwood,107032000
+Saint Kitts and Nevis,fuelwood,0
+Dominica,fuelwood,7517
+Antigua and Barbuda,fuelwood,0
+Aruba,fuelwood,1748
+Grenada,fuelwood,0
+Saint Vincent and the Grenadines,fuelwood,7444
+Saint Lucia,fuelwood,9850
+Barbados,fuelwood,5004
+Bahamas,fuelwood,32896
+Trinidad and Tobago,fuelwood,32822
+Jamaica,fuelwood,545051
+Puerto Rico,fuelwood,0
+Dominican Republic,fuelwood,912654
+Haiti,fuelwood,2041135
+Cuba,fuelwood,1141000
+Bermuda,fuelwood,0
+Belize,fuelwood,126000
+Panama,fuelwood,1128086
+Costa Rica,fuelwood,3377039
+Nicaragua,fuelwood,6096562
+El Salvador,fuelwood,4225392
+Honduras,fuelwood,8575298
+Guatemala,fuelwood,18058663
+Colombia,fuelwood,8826000
+Mexico,fuelwood,38828976
+Uruguay,fuelwood,2430000
+Paraguay,fuelwood,6575874
+Ecuador,fuelwood,4940180
+Peru,fuelwood,7338000
+Suriname,fuelwood,97100
+Guyana,fuelwood,848096
+Venezuela (Bolivarian Republic of),fuelwood,4054562
+Algeria,fuelwood,8176480
+Egypt,fuelwood,17511447
+Iran (Islamic Republic of),fuelwood,53000
+Iraq,fuelwood,118000
+Bahrain,fuelwood,6595
+Qatar,fuelwood,4976
+Oman,fuelwood,38976
+Kuwait,fuelwood,18406
+Lebanon,fuelwood,18866
+Jordan,fuelwood,302465
+Israel,fuelwood,2043
+United Arab Emirates,fuelwood,18002
+Yemen,fuelwood,440785
+Saudi Arabia,fuelwood,247428
+Morocco,fuelwood,6785689
+Malta,fuelwood,0
+Djibouti,fuelwood,355953
+Libya,fuelwood,951892
+Tunisia,fuelwood,2184607
+Canada,fuelwood,3211000
+United States of America,fuelwood,40436676
+Afghanistan,fuelwood,1641290
+Bangladesh,fuelwood,27286834
+Maldives,fuelwood,15398
+Sri Lanka,fuelwood,5211585
+India,fuelwood,309306678
+Bhutan,fuelwood,4845000
+Nepal,fuelwood,12525685
+Pakistan,fuelwood,29660000
+Sao Tome and Principe,fuelwood,107426
+Comoros,fuelwood,265913
+Equatorial Guinea,fuelwood,447000
+Congo,fuelwood,1335790
+Central African Republic,fuelwood,2000000
+Burundi,fuelwood,9259292
+Rwanda,fuelwood,5000000
+Chad,fuelwood,7070029
+Cameroon,fuelwood,9905983
+Democratic Republic of the Congo,fuelwood,76602030
+Eritrea,fuelwood,920987
+South Sudan,fuelwood,0
+Ethiopia,fuelwood,101274300
+Kenya,fuelwood,26400000
+Nigeria,fuelwood,63214728
+South Africa,fuelwood,12000000
+Swaziland,fuelwood,1314000
+Mauritius,fuelwood,9000
+Lesotho,fuelwood,2091851
+Botswana,fuelwood,683295
+Namibia,fuelwood,1327638
+Zambia,fuelwood,21823000
+Zimbabwe,fuelwood,8708870
+Malawi,fuelwood,5405340
+Madagascar,fuelwood,13100000
+Angola,fuelwood,4009338
+Mozambique,fuelwood,16724000
+Sudan,fuelwood,18775768
+Uganda,fuelwood,39636358
+United Republic of Tanzania,fuelwood,22835697
+Cabo Verde,fuelwood,186000
+Gabon,fuelwood,1070000
+Guinea-Bissau,fuelwood,2599954
+Gambia,fuelwood,694274
+Mauritania,fuelwood,1836252
+Liberia,fuelwood,7008341
+Sierra Leone,fuelwood,5581619
+Togo,fuelwood,4424000
+Benin,fuelwood,6274503
+Guinea,fuelwood,11958951
+Senegal,fuelwood,5427205
+Mali,fuelwood,5326348
+Burkina Faso,fuelwood,12784959
+Niger,fuelwood,9875857
+Cote d'Ivoire,fuelwood,8946768
+Ghana,fuelwood,37790966
+Somalia,fuelwood,13501108
diff --git a/data/wood_net_imports.csv b/data/wood_net_imports.csv
new file mode 100644
index 0000000000000000000000000000000000000000..d96a71c22dc0fc674496f005ddd8bfff0c7e144e
--- /dev/null
+++ b/data/wood_net_imports.csv
@@ -0,0 +1,387 @@
+country,type,harvest,net_import
+New Zealand,roundwood,24490000,-10742795
+Australia,roundwood,25561000,-1446019
+New Caledonia,roundwood,14700,589
+China,roundwood,161809500,35366442
+Tuvalu,roundwood,0,64
+Palau,roundwood,0,632
+Marshall Islands,roundwood,0,0
+Kiribati,roundwood,0,6
+Micronesia (Federated States of),roundwood,0,493
+Tonga,roundwood,2000,14586
+Samoa,roundwood,5700,617
+Vanuatu,roundwood,38000,-102
+Brunei Darussalam,roundwood,107300,-113
+Solomon Islands,roundwood,1624000,-1594000
+Macau,roundwood,0,0
+Fiji,roundwood,800000,-3714
+Lao People's Democratic Republic,roundwood,1032000,-206734
+Papua New Guinea,roundwood,4476000,-2784211
+French Polynesia,roundwood,1000,2466
+Democratic People's Republic of Korea,roundwood,1500000,-9768
+Indonesia,roundwood,54105500,25888
+Japan,roundwood,17193000,4691141
+Mongolia,roundwood,49000,-167
+Myanmar,roundwood,4171000,-1864473
+Singapore,roundwood,0,13648
+Malaysia,roundwood,19303000,-4456106
+Philippines,roundwood,3593000,6731
+Republic of Korea,roundwood,3466000,4216740
+Thailand,roundwood,14600000,259173
+Cambodia,roundwood,288000,-5638
+Viet Nam,roundwood,24200000,782678
+Albania,roundwood,80000,1393
+Armenia,roundwood,1000,-13
+Slovakia,roundwood,9089175,-1851697
+Serbia,roundwood,1413000,72310
+Bulgaria,roundwood,3011000,-447176
+Hungary,roundwood,2746292,-612081
+Czechia,roundwood,14771000,-764000
+Romania,roundwood,10548052,242426
+Montenegro,roundwood,208001,-13660
+Kosovo,roundwood,0,0
+Slovenia,roundwood,1841404,-372103
+The former Yugoslav Republic of Macedonia,roundwood,101000,33814
+Republic of Moldova,roundwood,43000,36452
+Bosnia and Herzegovina,roundwood,2461000,-105910
+Croatia,roundwood,3421000,-572000
+Luxembourg,roundwood,420593,691228
+Ireland,roundwood,2436975,-232948
+Belgium,roundwood,4113900,3337797
+Netherlands,roundwood,790593,-270800
+France,roundwood,29634241,-4974701
+Switzerland,roundwood,3439389,-507648
+Austria,roundwood,13281444,7085869
+Germany,roundwood,47136000,3930044
+Italy,roundwood,2647228,3151460
+Georgia,roundwood,105006,15794
+Turkmenistan,roundwood,0,0
+Kyrgyzstan,roundwood,9300,3300
+Tajikistan,roundwood,0,0
+Azerbaijan,roundwood,3300,10690
+Belarus,roundwood,8072600,-2182394
+Kazakhstan,roundwood,73000,85136
+Uzbekistan,roundwood,8000,270257
+Poland,roundwood,31343001,704011
+Estonia,roundwood,4841309,-1949726
+Latvia,roundwood,10221818,-3721174
+Lithuania,roundwood,5153862,-1048474
+Russian Federation,roundwood,161595000,-20981640
+Greenland,roundwood,0,72
+Iceland,roundwood,0,150
+Norway,roundwood,8322430,423279
+Finland,roundwood,45419600,5772749
+Denmark,roundwood,1695300,151918
+Sweden,roundwood,66300000,5059465
+Portugal,roundwood,9048360,-145539
+Spain,roundwood,10969399,507246
+Cyprus,roundwood,5330,357
+Greece,roundwood,336481,131693
+Turkey,roundwood,15695000,1166600
+Ukraine,roundwood,7536000,-2913700
+United Kingdom,roundwood,8337163,-121323
+Argentina,roundwood,10635024,-45742
+Bolivia (Plurinational State of),roundwood,943000,-5909
+Chile,roundwood,34560000,-34835
+Brazil,roundwood,128400000,-76
+Saint Kitts and Nevis,roundwood,0,351
+Dominica,roundwood,0,831
+Antigua and Barbuda,roundwood,0,377
+Aruba,roundwood,0,2667
+Grenada,roundwood,0,0
+Saint Vincent and the Grenadines,roundwood,0,5507
+Saint Lucia,roundwood,0,925
+Barbados,roundwood,6000,1158
+Bahamas,roundwood,17000,7096
+Trinidad and Tobago,roundwood,65000,170
+Jamaica,roundwood,154805,12134
+Puerto Rico,roundwood,0,0
+Dominican Republic,roundwood,9700,12587
+Haiti,roundwood,239000,24186
+Cuba,roundwood,589500,16142
+Bermuda,roundwood,0,0
+Belize,roundwood,41000,5542
+Panama,roundwood,176000,-43651
+Costa Rica,roundwood,1326000,-113515
+Nicaragua,roundwood,108000,-3293
+El Salvador,roundwood,682000,-23146
+Honduras,roundwood,620000,-1858
+Guatemala,roundwood,692000,-9547
+Colombia,roundwood,3550000,-17575
+Mexico,roundwood,5041000,13930
+Uruguay,roundwood,9402000,-466000
+Paraguay,roundwood,4044000,-16054
+Ecuador,roundwood,2091000,-110442
+Peru,roundwood,1352000,355
+Suriname,roundwood,246300,-48470
+Guyana,roundwood,516000,-115574
+Venezuela (Bolivarian Republic of),roundwood,1850000,401
+Algeria,roundwood,144600,10086
+Egypt,roundwood,268000,131388
+Iran (Islamic Republic of),roundwood,697000,16945
+Iraq,roundwood,59000,1353
+Bahrain,roundwood,0,3995
+Qatar,roundwood,0,3543
+Oman,roundwood,0,12331
+Kuwait,roundwood,0,107
+Lebanon,roundwood,18000,7986
+Jordan,roundwood,4000,1389
+Israel,roundwood,24957,4696
+United Arab Emirates,roundwood,0,158934
+Yemen,roundwood,0,1271
+Saudi Arabia,roundwood,0,23430
+Morocco,roundwood,372100,238170
+Malta,roundwood,0,296
+Djibouti,roundwood,0,1416
+Libya,roundwood,116000,16306
+Tunisia,roundwood,218000,14586
+Canada,roundwood,138802000,726000
+United States of America,roundwood,336135000,-10766210
+Afghanistan,roundwood,1760000,-339
+Bangladesh,roundwood,420000,16733
+Maldives,roundwood,0,44
+Sri Lanka,roundwood,698000,594
+India,roundwood,48759000,5296057
+Bhutan,roundwood,190000,-654
+Nepal,roundwood,1300000,-8
+Pakistan,roundwood,4060000,141869
+Sao Tome and Principe,roundwood,27200,-80
+Comoros,roundwood,24650,0
+Equatorial Guinea,roundwood,309849,-226811
+Congo,roundwood,2035524,-606082
+Central African Republic,roundwood,632283,-147893
+Burundi,roundwood,883000,-2493
+Rwanda,roundwood,1211927,1158
+Chad,roundwood,761000,508
+Cameroon,roundwood,2348000,-701454
+Democratic Republic of the Congo,roundwood,4531539,-162422
+Eritrea,roundwood,5654,104
+South Sudan,roundwood,543250,-2368
+Ethiopia,roundwood,2935000,1000
+Kenya,roundwood,1133000,14092
+Nigeria,roundwood,9418000,-109896
+South Africa,roundwood,16988569,-169917
+Swaziland,roundwood,934000,1117
+Mauritius,roundwood,5000,27051
+Lesotho,roundwood,0,2995
+Botswana,roundwood,105000,8360
+Namibia,roundwood,0,4266
+Zambia,roundwood,1455000,-3220
+Zimbabwe,roundwood,518400,34
+Malawi,roundwood,1480000,-5648
+Madagascar,roundwood,280911,3091
+Angola,roundwood,1127000,-2684
+Mozambique,roundwood,1416000,-235956
+Sudan,roundwood,1629750,-7102
+Uganda,roundwood,4093000,645
+United Republic of Tanzania,roundwood,2314000,-8995
+Cabo Verde,roundwood,0,1274
+Gabon,roundwood,1800000,-905000
+Guinea-Bissau,roundwood,131910,-14963
+Gambia,roundwood,131700,-18350
+Mauritania,roundwood,32000,18058
+Liberia,roundwood,480000,-7729
+Sierra Leone,roundwood,123600,-14885
+Togo,roundwood,166000,-117421
+Benin,roundwood,467000,-86400
+Guinea,roundwood,651000,-149183
+Senegal,roundwood,788000,3761
+Mali,roundwood,412900,-351
+Burkina Faso,roundwood,1171000,-537
+Niger,roundwood,701000,543
+Cote d'Ivoire,roundwood,2356000,-145146
+Ghana,roundwood,1950000,-141277
+Somalia,roundwood,110000,1187
+New Zealand,fuelwood,0,-96
+Australia,fuelwood,4853000,-260
+New Caledonia,fuelwood,11818,1
+China,fuelwood,188823052,2840
+Tuvalu,fuelwood,0,0
+Palau,fuelwood,0,0
+Marshall Islands,fuelwood,0,0
+Kiribati,fuelwood,2937,121
+Micronesia (Federated States of),fuelwood,2444,215
+Tonga,fuelwood,2134,3
+Samoa,fuelwood,70000,14
+Vanuatu,fuelwood,91000,-596
+Brunei Darussalam,fuelwood,11679,68
+Solomon Islands,fuelwood,128451,0
+Macau,fuelwood,0,0
+Fiji,fuelwood,37000,-11
+Lao People's Democratic Republic,fuelwood,5947951,-189
+Papua New Guinea,fuelwood,5533000,0
+French Polynesia,fuelwood,4417,96
+Democratic People's Republic of Korea,fuelwood,5986734,8
+Indonesia,fuelwood,59743242,53127
+Japan,fuelwood,87729,-162
+Mongolia,fuelwood,763361,42
+Myanmar,fuelwood,38286000,-19
+Singapore,fuelwood,32799,-186
+Malaysia,fuelwood,2810489,-4870
+Philippines,fuelwood,12361935,0
+Republic of Korea,fuelwood,518000,47
+Thailand,fuelwood,19301440,6103
+Cambodia,fuelwood,8442343,0
+Viet Nam,fuelwood,21500000,-32741
+Albania,fuelwood,350000,-56300
+Armenia,fuelwood,1750000,0
+Slovakia,fuelwood,509893,-61769
+Serbia,fuelwood,6223000,-9600
+Bulgaria,fuelwood,2657000,-192375
+Hungary,fuelwood,2993983,2414
+Czechia,fuelwood,1965000,-26000
+Romania,fuelwood,2563588,-60632
+Montenegro,fuelwood,707000,-2913
+Kosovo,fuelwood,0,0
+Slovenia,fuelwood,1104045,-165052
+The former Yugoslav Republic of Macedonia,fuelwood,530000,42935
+Republic of Moldova,fuelwood,308800,2364
+Bosnia and Herzegovina,fuelwood,1418000,-488170
+Croatia,fuelwood,1056000,-240000
+Luxembourg,fuelwood,81576,-2917
+Ireland,fuelwood,181021,4266
+Belgium,fuelwood,713525,44488
+Netherlands,fuelwood,290000,-2000
+France,fuelwood,26173567,-770532
+Switzerland,fuelwood,1498647,-13625
+Austria,fuelwood,4549512,535257
+Germany,fuelwood,27296000,282693
+Italy,fuelwood,11429000,951000
+Georgia,fuelwood,733000,0
+Turkmenistan,fuelwood,10000,0
+Kyrgyzstan,fuelwood,36600,-151
+Tajikistan,fuelwood,90000,0
+Azerbaijan,fuelwood,3200,0
+Belarus,fuelwood,2291600,-4119
+Kazakhstan,fuelwood,272000,70
+Uzbekistan,fuelwood,22000,0
+Poland,fuelwood,4124416,-113720
+Estonia,fuelwood,2195874,-167019
+Latvia,fuelwood,2312000,-1327464
+Lithuania,fuelwood,1942998,-61339
+Russian Federation,fuelwood,13904000,-193293
+Greenland,fuelwood,0,9
+Iceland,fuelwood,0,34
+Norway,fuelwood,2120652,171664
+Finland,fuelwood,6705000,133296
+Denmark,fuelwood,1362200,136461
+Sweden,fuelwood,5900000,475717
+Portugal,fuelwood,600000,-1162
+Spain,fuelwood,5120000,-57310
+Cyprus,fuelwood,3628,237
+Greece,fuelwood,711481,296177
+Turkey,fuelwood,4902000,241988
+Ukraine,fuelwood,8609600,-737859
+United Kingdom,fuelwood,1381100,-155755
+Argentina,fuelwood,4374522,-161
+Bolivia (Plurinational State of),fuelwood,2349841,0
+Chile,fuelwood,12655000,0
+Brazil,fuelwood,107032000,114
+Saint Kitts and Nevis,fuelwood,0,0
+Dominica,fuelwood,7517,-24
+Antigua and Barbuda,fuelwood,0,0
+Aruba,fuelwood,1748,1
+Grenada,fuelwood,0,0
+Saint Vincent and the Grenadines,fuelwood,7444,0
+Saint Lucia,fuelwood,9850,0
+Barbados,fuelwood,5004,-22
+Bahamas,fuelwood,32896,11
+Trinidad and Tobago,fuelwood,32822,53
+Jamaica,fuelwood,545051,2
+Puerto Rico,fuelwood,0,0
+Dominican Republic,fuelwood,912654,-7
+Haiti,fuelwood,2041135,7
+Cuba,fuelwood,1141000,0
+Bermuda,fuelwood,0,0
+Belize,fuelwood,126000,79
+Panama,fuelwood,1128086,0
+Costa Rica,fuelwood,3377039,-24
+Nicaragua,fuelwood,6096562,-783
+El Salvador,fuelwood,4225392,1
+Honduras,fuelwood,8575298,-407
+Guatemala,fuelwood,18058663,0
+Colombia,fuelwood,8826000,28
+Mexico,fuelwood,38828976,-8366
+Uruguay,fuelwood,2430000,0
+Paraguay,fuelwood,6575874,-57
+Ecuador,fuelwood,4940180,0
+Peru,fuelwood,7338000,0
+Suriname,fuelwood,97100,-396
+Guyana,fuelwood,848096,-43
+Venezuela (Bolivarian Republic of),fuelwood,4054562,0
+Algeria,fuelwood,8176480,2
+Egypt,fuelwood,17511447,530
+Iran (Islamic Republic of),fuelwood,53000,972
+Iraq,fuelwood,118000,0
+Bahrain,fuelwood,6595,249
+Qatar,fuelwood,4976,372
+Oman,fuelwood,38976,-33
+Kuwait,fuelwood,18406,98
+Lebanon,fuelwood,18866,1594
+Jordan,fuelwood,302465,144
+Israel,fuelwood,2043,87
+United Arab Emirates,fuelwood,18002,4622
+Yemen,fuelwood,440785,91
+Saudi Arabia,fuelwood,247428,1611
+Morocco,fuelwood,6785689,169
+Malta,fuelwood,0,25
+Djibouti,fuelwood,355953,23
+Libya,fuelwood,951892,0
+Tunisia,fuelwood,2184607,24
+Canada,fuelwood,3211000,44879
+United States of America,fuelwood,40436676,-173728
+Afghanistan,fuelwood,1641290,13
+Bangladesh,fuelwood,27286834,164
+Maldives,fuelwood,15398,179
+Sri Lanka,fuelwood,5211585,-157
+India,fuelwood,309306678,6971
+Bhutan,fuelwood,4845000,0
+Nepal,fuelwood,12525685,13
+Pakistan,fuelwood,29660000,-14
+Sao Tome and Principe,fuelwood,107426,2
+Comoros,fuelwood,265913,0
+Equatorial Guinea,fuelwood,447000,0
+Congo,fuelwood,1335790,0
+Central African Republic,fuelwood,2000000,0
+Burundi,fuelwood,9259292,-11
+Rwanda,fuelwood,5000000,-20
+Chad,fuelwood,7070029,5
+Cameroon,fuelwood,9905983,-84
+Democratic Republic of the Congo,fuelwood,76602030,-16
+Eritrea,fuelwood,920987,0
+South Sudan,fuelwood,18775768,271
+Ethiopia,fuelwood,101274300,210
+Kenya,fuelwood,26400000,2590
+Nigeria,fuelwood,63214728,406
+South Africa,fuelwood,12000000,268000
+Swaziland,fuelwood,1314000,-251529
+Mauritius,fuelwood,9000,-3
+Lesotho,fuelwood,2091851,772
+Botswana,fuelwood,683295,130
+Namibia,fuelwood,1327638,-19964
+Zambia,fuelwood,21823000,-78
+Zimbabwe,fuelwood,8708870,65
+Malawi,fuelwood,5405340,-1002
+Madagascar,fuelwood,13100000,0
+Angola,fuelwood,4009338,5
+Mozambique,fuelwood,16724000,246
+Sudan,fuelwood,18775768,271
+Uganda,fuelwood,39636358,-13
+United Republic of Tanzania,fuelwood,22835697,-2670
+Cabo Verde,fuelwood,186000,55
+Gabon,fuelwood,1070000,-2
+Guinea-Bissau,fuelwood,2599954,0
+Gambia,fuelwood,694274,4
+Mauritania,fuelwood,1836252,0
+Liberia,fuelwood,7008341,0
+Sierra Leone,fuelwood,5581619,-32
+Togo,fuelwood,4424000,-32
+Benin,fuelwood,6274503,1
+Guinea,fuelwood,11958951,0
+Senegal,fuelwood,5427205,2
+Mali,fuelwood,5326348,-505
+Burkina Faso,fuelwood,12784959,6
+Niger,fuelwood,9875857,0
+Cote d'Ivoire,fuelwood,8946768,-1238
+Ghana,fuelwood,37790966,-8027
+Somalia,fuelwood,13501108,-245
diff --git a/debug_config.properties b/debug_config.properties
index 1432cdd4b78a712a13d35ca2707c5ba6f0039848..9b72a2af42b62229ddbf459009c81797d6f4e6e3 100644
--- a/debug_config.properties
+++ b/debug_config.properties
@@ -7,10 +7,10 @@ OUTPUT_DIR = C:/Users/Bart/git/plumv2/output
BASE_YEAR=2010
START_TIMESTEP=0
TIMESTEP_SIZE=1
-END_TIMESTEP=90
+END_TIMESTEP=10
IS_CALIBRATION_RUN=true
-END_FIRST_STAGE_CALIBRATION=2
+END_FIRST_STAGE_CALIBRATION=0
GENERATE_NEW_YIELD_CLUSTERS=false
@@ -20,12 +20,14 @@ DEBUG_LIMIT_COUNTRIES=true
DEBUG_COUNTRY_NAME=Brazil
GAMS_COUNTRY_TO_SAVE=Brazil
-INIT_WOOD_PRICE=0.2
+INIT_WOOD_PRICE=0.0
INIT_WOOD_STOCK=1000
INIT_CARBON_PRICE=0.0
FOREST_ROTATION_PERIOD=30
INFRASTRUCTURE_EXPANSION_COST=0.0
-MANAGED_FOREST_INCREASE_COST=0.03
-MANAGED_FOREST_DECREASE_COST=0.03
-FOREST_MANAGEMENT_COST=0.05
-FOREST_MAX_CHANGE=0.02
\ No newline at end of file
+MANAGED_FOREST_INCREASE_COST=0.2
+MANAGED_FOREST_DECREASE_COST=0.2
+FOREST_MANAGEMENT_COST=0.2
+FOREST_MAX_CHANGE=0.02
+LAND_CHANGE_COST=0.7
+NATURAL_AREA_VALUE=0.01
\ No newline at end of file
diff --git a/src/ac/ed/lurg/InternationalMarket.java b/src/ac/ed/lurg/InternationalMarket.java
index bd5735fae189417f01da2323e43dfaad28e017e9..06deaee9eedbf6058455276829a4d8989ae7db2d 100644
--- a/src/ac/ed/lurg/InternationalMarket.java
+++ b/src/ac/ed/lurg/InternationalMarket.java
@@ -17,9 +17,11 @@ import ac.ed.lurg.country.AbstractCountryAgent;
import ac.ed.lurg.country.GlobalPrice;
import ac.ed.lurg.country.StockReader;
import ac.ed.lurg.landuse.CropUsageData;
+import ac.ed.lurg.landuse.WoodUsageData;
import ac.ed.lurg.shock.PriceShockManager;
import ac.ed.lurg.types.CropToDoubleMap;
import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.WoodType;
import ac.ed.lurg.utils.LogWriter;
public class InternationalMarket {
@@ -137,12 +139,15 @@ public class InternationalMarket {
double totalWoodExport = 0;
double totalWoodProduction = 0;
for (AbstractCountryAgent ca : countryAgents) {
- totalWoodProduction += ca.getWoodUsageData().getHarvest();
- double netImport = ca.getWoodUsageData().getNetImport();
- if (netImport >= 0) {
- totalWoodImport += netImport;
- } else {
- totalWoodExport += -netImport;
+ Map<WoodType, WoodUsageData> woodUsageMap = ca.getWoodUsageData();
+ for (WoodUsageData woodUsage : woodUsageMap.values()) {
+ totalWoodProduction += woodUsage.getHarvest();
+ double netImport = woodUsage.getNetImport();
+ if (netImport >= 0) {
+ totalWoodImport += netImport;
+ } else {
+ totalWoodExport += -netImport;
+ }
}
}
totalWoodProduction = Math.max(totalWoodProduction, 0.0000001); // avoid division by 0
diff --git a/src/ac/ed/lurg/ModelConfig.java b/src/ac/ed/lurg/ModelConfig.java
index f1d0c3031e3818a984546c5729f688fa1edfe313..f62a7b207ccd4736ebc2051d44baef4d46c1eafd 100755
--- a/src/ac/ed/lurg/ModelConfig.java
+++ b/src/ac/ed/lurg/ModelConfig.java
@@ -250,6 +250,9 @@ public class ModelConfig {
public static final String FOREST_DIR = SPATIAL_DATA_DIR + File.separator + "forestry";
//public static final String CARBON_FLUX_FILE = FOREST_DIR + File.separator + "carbon_flux_";
public static final String WOOD_YIELD_FILENAME = "wood_yield.out";
+ public static final String TIMBER_FOREST_GROWTH_FILENAME = "x1.out";
+ public static final String CARBON_FOREST_GROWTH_FILENAME = "x2.out";
+ public static final String NATURAL_FOREST_GROWTH_FILENAME = "x3.out";
// Output
public static final String LAND_COVER_OUTPUT_FILE = OUTPUT_DIR + File.separator + "lc.txt";
@@ -275,7 +278,7 @@ public class ModelConfig {
public static final boolean CHANGE_YIELD_DATA_YEAR = IS_CALIBRATION_RUN ? false : getBooleanProperty("CHANGE_YIELD_DATA_YEAR", true);
public static final String CLUSTERED_YIELD_FILE = getProperty("CLUSTERED_YIELD_FILE", CALIB_DIR + File.separator + "cluster.asc");
public static final boolean GENERATE_NEW_YIELD_CLUSTERS = getBooleanProperty("GENERATE_NEW_YIELD_CLUSTERS", IS_CALIBRATION_RUN);
- public static final String SERIALIZED_WOOD_AND_CARBON_PROD_FILE = CALIB_DIR + File.separator + "countryWoodAndCarbonProd.ser";
+ public static final String SERIALIZED_WOOD_USAGE_FILE = CALIB_DIR + File.separator + "countryWoodUsage.ser";
// Temporal configuration
public static final int START_TIMESTEP = getIntProperty("START_TIMESTEP", 0);
@@ -440,4 +443,10 @@ public class ModelConfig {
public static final double FOREST_MAX_CHANGE = getDoubleProperty("FOREST_MAX_CHANGE", 0.02);
public static final double ROUNDWOOD_DEMAND_ELASTICITY = getDoubleProperty("ROUNDWOOD_DEMAND_ELASTICITY", 0.5);
public static final double FUELWOOD_DEMAND_ELASTICITY = getDoubleProperty("FUELDWOOD_DEMAND_ELASTICITY", -0.5);
+ public static final double FOREST_ESTABLISHMENT_COST = getDoubleProperty("FOREST_ESTABLISHMENT_COST", 0.2);
+ public static final double NATURAL_AREA_VALUE = getDoubleProperty("NATURAL_AREA_VALUE", 0.1); // $1000/ha
+ public static final String WOOD_NET_IMPORTS_FILENAME = getProperty("WOOD_NET_IMPORTS_FILENAME", "wood_net_imports.csv");
+ public static final String WOOD_NET_IMPORTS_FILE = getProperty("WOOD_NET_IMPORTS_FILE", DATA_DIR + File.separator + WOOD_NET_IMPORTS_FILENAME);
+ public static final double WOOD_BIOMASS_CONVERSION_FACTOR = getDoubleProperty("WOOD_BIOMASS_CONVERSION_FACTOR", 2.688e-7); // m3 to MtC-eq p.16 [https://doi.org/10.5194/gmd-13-5425-2020]
+ public static final int FOREST_INIT_AGE = getIntProperty("FOREST_INIT_AGE", 60);
}
diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java
index 357c675b64f4841486f5938556f2e25fd37ce19f..487cfa3846a1c9c00dc3fecdb08504727f97361f 100644
--- a/src/ac/ed/lurg/ModelMain.java
+++ b/src/ac/ed/lurg/ModelMain.java
@@ -27,6 +27,10 @@ import ac.ed.lurg.demand.CalorieManager;
import ac.ed.lurg.demand.DemandManagerFromFile;
import ac.ed.lurg.demand.DemandManagerSSP;
import ac.ed.lurg.demand.ElasticDemandManager;
+import ac.ed.lurg.forestry.CarbonForestGrowthDataReader;
+import ac.ed.lurg.forestry.ForestGrowthRasterSet;
+import ac.ed.lurg.forestry.NaturalForestGrowthDataReader;
+import ac.ed.lurg.forestry.TimberForestGrowthDataReader;
import ac.ed.lurg.landuse.CarbonFluxRasterSet;
import ac.ed.lurg.landuse.CarbonFluxReader;
import ac.ed.lurg.landuse.ConversionCostReader;
@@ -44,6 +48,8 @@ import ac.ed.lurg.landuse.LandUseItem;
import ac.ed.lurg.landuse.MaxCropAreaReader;
import ac.ed.lurg.landuse.ProtectedAreasReader;
import ac.ed.lurg.landuse.RunOffReader;
+import ac.ed.lurg.landuse.WoodUsageData;
+import ac.ed.lurg.landuse.WoodUsageReader;
import ac.ed.lurg.landuse.WoodYieldRasterSet;
import ac.ed.lurg.landuse.WoodYieldReader;
import ac.ed.lurg.output.LandUseOutputer;
@@ -51,6 +57,7 @@ import ac.ed.lurg.output.LpjgOutputer;
import ac.ed.lurg.types.CommodityType;
import ac.ed.lurg.types.CropType;
import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.types.WoodType;
import ac.ed.lurg.utils.DoubleMap;
import ac.ed.lurg.utils.FileWriterHelper;
import ac.ed.lurg.utils.LogWriter;
@@ -469,7 +476,7 @@ public class ModelMain {
if (ModelConfig.IS_CALIBRATION_RUN) {
serializeLandUse(landUseRaster);
countryAgents.serializeCropUsageForAll();
- countryAgents.serializeTimberAndCarbonProdForAll();
+ countryAgents.serializeWoodUsageForAll();
internationalMarket.serializeGlobalPrices();
}
@@ -537,7 +544,7 @@ public class ModelMain {
public void createCountryAgents(Collection<CompositeCountry> countryGrouping) {
countryAgents = new CountryAgentManager(compositeCountryManager, demandManager, countryBoundaryRaster, internationalMarket, clusterIdRaster, globalLandUseRaster);
Map<CompositeCountry, Map<CropType, CropUsageData>> cropUsageDataMap = getInitialCropUsageData();
- Map<CompositeCountry, Double[]> woodAndCarbonProdMap = getInitialWoodAndCarbonProd();
+ Map<CompositeCountry, Map<WoodType, WoodUsageData>> woodAndCarbonProdMap = getInitialWoodUsageData();
RasterSet<LandUseItem> initLU = getInitialLandUse();
for (CompositeCountry cc : countryGrouping) {
@@ -566,33 +573,30 @@ public class ModelMain {
return cropUsageDataMap;
}
- private Map<CompositeCountry, Double[]> getInitialWoodAndCarbonProd() {
- Map<CompositeCountry, Double[]> initWoodAndCarbonProd;
+ private Map<CompositeCountry, Map<WoodType, WoodUsageData>> getInitialWoodUsageData() {
+ Map<CompositeCountry, Map<WoodType, WoodUsageData>> woodUsageDataMap;
if (ModelConfig.IS_CALIBRATION_RUN) {
- initWoodAndCarbonProd = new HashMap<CompositeCountry, Double[]>();
- Double[] initProdValues = {0.0, 0.0}; // TODO Read in
- for (CompositeCountry cc : compositeCountryManager.getAll()) {
- initWoodAndCarbonProd.put(cc, initProdValues);
- }
+ woodUsageDataMap = new WoodUsageReader(compositeCountryManager).getWoodUsageData();
+
} else {
- initWoodAndCarbonProd = deserialiseWoodAndCarbonProd();
+ woodUsageDataMap = deserializeWoodUsage();
}
- return initWoodAndCarbonProd;
+ return woodUsageDataMap;
}
@SuppressWarnings("unchecked")
- private Map<CompositeCountry, Double[]> deserialiseWoodAndCarbonProd() {
+ private Map<CompositeCountry, Map<WoodType, WoodUsageData>> deserializeWoodUsage() {
try {
- Map<CompositeCountry, Double[]> initTimberAndCarbonProd;
- FileInputStream fileIn = new FileInputStream(ModelConfig.SERIALIZED_WOOD_AND_CARBON_PROD_FILE);
+ Map<CompositeCountry, Map<WoodType, WoodUsageData>> woodUsageDataMap;
+ FileInputStream fileIn = new FileInputStream(ModelConfig.SERIALIZED_WOOD_USAGE_FILE);
ObjectInputStream in = new ObjectInputStream(fileIn);
- initTimberAndCarbonProd = (Map<CompositeCountry, Double[]>) in.readObject();
+ woodUsageDataMap = (Map<CompositeCountry, Map<WoodType, WoodUsageData>>) in.readObject();
in.close();
fileIn.close();
- LogWriter.println("Deserialized " + ModelConfig.SERIALIZED_WOOD_AND_CARBON_PROD_FILE);
- return initTimberAndCarbonProd;
+ LogWriter.println("Deserialized " + ModelConfig.SERIALIZED_WOOD_USAGE_FILE);
+ return woodUsageDataMap;
} catch (IOException i) {
- LogWriter.printlnError("Problem deserializing " + ModelConfig.SERIALIZED_WOOD_AND_CARBON_PROD_FILE);
+ LogWriter.printlnError("Problem deserializing " + ModelConfig.SERIALIZED_WOOD_USAGE_FILE);
LogWriter.print(i);
return null;
} catch (ClassNotFoundException c) {
@@ -715,6 +719,15 @@ public class ModelMain {
new WoodYieldReader(woodYieldData).getRasterDataFromFile(timestep.getYearSubDir(ModelConfig.FOREST_DIR) + File.separator + ModelConfig.WOOD_YIELD_FILENAME);
return woodYieldData;
}
+
+ private ForestGrowthRasterSet getForestGrowthData(Timestep timestep) {
+ ForestGrowthRasterSet forestGrowthData = new ForestGrowthRasterSet(desiredProjection);
+ LogWriter.println("Reading forest growth data");
+ new TimberForestGrowthDataReader(forestGrowthData).getRasterDataFromFile(timestep.getYearSubDir(ModelConfig.FOREST_DIR) + File.separator + ModelConfig.TIMBER_FOREST_GROWTH_FILENAME);
+ new CarbonForestGrowthDataReader(forestGrowthData).getRasterDataFromFile(timestep.getYearSubDir(ModelConfig.FOREST_DIR) + File.separator + ModelConfig.CARBON_FOREST_GROWTH_FILENAME);
+ new NaturalForestGrowthDataReader(forestGrowthData).getRasterDataFromFile(timestep.getYearSubDir(ModelConfig.FOREST_DIR) + File.separator + ModelConfig.NATURAL_FOREST_GROWTH_FILENAME);
+ return forestGrowthData;
+ }
/** Ugly in situ update of currentIrrigationData, better if IrrigationRasterSets were handled more immutably */
private void getUpdateIrrigationData(Timestep timestep, YieldRaster yieldSurfaces) {
diff --git a/src/ac/ed/lurg/country/AbstractCountryAgent.java b/src/ac/ed/lurg/country/AbstractCountryAgent.java
index 58d7adf9be0182676998634ca4bc135db80d7243..df4f033a96f633126f273b336e5714390f14bc14 100644
--- a/src/ac/ed/lurg/country/AbstractCountryAgent.java
+++ b/src/ac/ed/lurg/country/AbstractCountryAgent.java
@@ -180,6 +180,6 @@ public abstract class AbstractCountryAgent {
abstract public double getNetCarbonFlux();
- abstract public WoodUsageData getWoodUsageData();
+ abstract public Map<WoodType, WoodUsageData> getWoodUsageData();
}
\ 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 184e433c2b8db64071c7847fdcd2c507b4f4e764..96a982ad70d4a74dd33c6fbd061dac6ab1ca51d4 100644
--- a/src/ac/ed/lurg/country/CountryAgent.java
+++ b/src/ac/ed/lurg/country/CountryAgent.java
@@ -25,6 +25,7 @@ import ac.ed.lurg.landuse.WoodYieldItem;
import ac.ed.lurg.types.CommodityType;
import ac.ed.lurg.types.CropType;
import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.types.WoodType;
import ac.ed.lurg.utils.DoubleMap;
import ac.ed.lurg.utils.LogWriter;
import ac.ed.lurg.utils.cluster.Cluster;
@@ -45,13 +46,13 @@ public class CountryAgent extends AbstractCountryAgent {
public CountryAgent(AbstractDemandManager demandManager,CompositeCountry country, RasterSet<LandUseItem> cropAreaRaster,
Map<CropType, CropUsageData> cropUsageData, Map<CropType, Double> tradeBarriers, RasterSet<IntegerRasterItem> yieldClusters,
- Map<CropType, Double> subsidyRates, Double[] woodAndCarbonProd) {
+ Map<CropType, Double> subsidyRates, Map<WoodType, WoodUsageData> countryWoodData) {
super(demandManager, country, tradeBarriers);
this.yieldClusters = yieldClusters;
this.subsidyRates = subsidyRates;
- GamsRasterOutput initialData = new GamsRasterOutput(cropAreaRaster, cropUsageData, woodAndCarbonProd);
+ GamsRasterOutput initialData = new GamsRasterOutput(cropAreaRaster, cropUsageData, countryWoodData);
previousGamsRasterOutput = initialData;
saveGamsGdxFiles = (ModelConfig.GAMS_COUNTRY_TO_SAVE != null && country.getName().equals(ModelConfig.GAMS_COUNTRY_TO_SAVE));
@@ -225,24 +226,23 @@ public class CountryAgent extends AbstractCountryAgent {
}
// Timber import/export constraints
- TradeConstraint timberTradeConstraint;
- {
- double baseTrade;
- if (ModelConfig.IS_CALIBRATION_RUN && currentTimestep.isInitialTimestep()) {
- baseTrade = -1;
- } else {
- baseTrade = getWoodUsageData().getNetImport();
- }
+ Map<WoodType, TradeConstraint> woodTradeConstraints = new HashMap<WoodType, TradeConstraint>();
+
+ for (Map.Entry<WoodType, WoodUsageData> entry : previousGamsRasterOutput.getWoodUsageData().entrySet()) {
+ WoodUsageData woodUsage = entry.getValue();
+ WoodType woodType = entry.getKey();
+ double baseTrade = woodUsage.getNetImport();
double changeUp = 0.0;
double changeDown = 0.0;
if (allowedImportChange > 0.0) {
changeUp = changeDown = baseTrade * allowedImportChange;
}
- timberTradeConstraint = new TradeConstraint(baseTrade + changeDown, baseTrade - changeUp);
+ woodTradeConstraints.put(woodType, new TradeConstraint(baseTrade + changeDown, baseTrade - changeUp));
}
GamsCountryInput countryLevelInputs = new GamsCountryInput(country, currentProjectedDemand, currentCountryPrices, importConstraints,
- previousGamsRasterOutput.getCropUsageData(), currentMinDemandFract, subsidyRates, currentCarbonPrice, carbonTradeConstraint, currentWoodDemand, currentTimberPrice, timberTradeConstraint);
+ previousGamsRasterOutput.getCropUsageData(), currentMinDemandFract, subsidyRates, currentCarbonPrice, carbonTradeConstraint, currentWoodDemand,
+ currentTimberPrice, woodTradeConstraints, previousGamsRasterOutput.getWoodUsageData());
GamsRasterInput input = new GamsRasterInput(currentTimestep, countryYieldSurfaces, previousGamsRasterOutput.getLandUses(), irrigData,
carbonFluxData, woodYieldData, conversionCosts, countryLevelInputs);
@@ -336,7 +336,7 @@ public class CountryAgent extends AbstractCountryAgent {
return previousGamsRasterOutput.getNetCarbonFlux();
}
- public WoodUsageData getWoodUsageData() {
+ public Map<WoodType, WoodUsageData> getWoodUsageData() {
return previousGamsRasterOutput.getWoodUsageData();
}
}
\ No newline at end of file
diff --git a/src/ac/ed/lurg/country/CountryAgentManager.java b/src/ac/ed/lurg/country/CountryAgentManager.java
index 56d07768a0a823926876dfa19c0fb1f9dd5744b4..1683b9ad1944f3b655d4ec30a493305425180ae4 100644
--- a/src/ac/ed/lurg/country/CountryAgentManager.java
+++ b/src/ac/ed/lurg/country/CountryAgentManager.java
@@ -21,10 +21,12 @@ import ac.ed.lurg.landuse.CropUsageData;
import ac.ed.lurg.landuse.IrrigationItem;
import ac.ed.lurg.landuse.IrrigationRasterSet;
import ac.ed.lurg.landuse.LandUseItem;
+import ac.ed.lurg.landuse.WoodUsageData;
import ac.ed.lurg.landuse.WoodYieldItem;
import ac.ed.lurg.landuse.WoodYieldRasterSet;
import ac.ed.lurg.types.CropType;
import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.types.WoodType;
import ac.ed.lurg.utils.DoubleMap;
import ac.ed.lurg.utils.LogWriter;
import ac.ed.lurg.yield.YieldRaster;
@@ -60,7 +62,7 @@ public class CountryAgentManager {
}
public void addForCountry(CompositeCountry cc, Map<CompositeCountry, Map<CropType, CropUsageData>> cropUsageDataMap,
- RasterSet<LandUseItem> initLU, Map<CompositeCountry, Double[]> woodAndCarbonProdMap) {
+ RasterSet<LandUseItem> initLU, Map<CompositeCountry, Map<WoodType, WoodUsageData>> woodUsageDataMap) {
if (ModelConfig.DEBUG_LIMIT_COUNTRIES) { // DEBUG code
if (!(cc.getName().equals(ModelConfig.DEBUG_COUNTRY_NAME)))
return;
@@ -78,7 +80,7 @@ public class CountryAgentManager {
LogWriter.println("Creating GAMS agent for: " + cc);
List<RasterKey> keys = countryBoundaryRaster.getKeysFor(cc);
Map<CropType, CropUsageData> countryCommodityData = cropUsageDataMap.get(cc);
- Double[] woodAndCarbonProd = woodAndCarbonProdMap.get(cc);
+ Map<WoodType, WoodUsageData> countryWoodData = woodUsageDataMap.get(cc);
if (countryCommodityData == null) {
LogWriter.printlnError("No commodities data for " + cc + ", so skipping");
@@ -92,7 +94,7 @@ public class CountryAgentManager {
else {
Map<CropType, Double> subsidyRates = subsidyRateManager.getSubsidyRates(cc);
- CountryAgent ca = new CountryAgent(demandManager, cc, initCountryLandUse, countryCommodityData, tradeBarriers, yieldClusters, subsidyRates, woodAndCarbonProd);
+ CountryAgent ca = new CountryAgent(demandManager, cc, initCountryLandUse, countryCommodityData, tradeBarriers, yieldClusters, subsidyRates, countryWoodData);
gamsCountryAgents.add(ca);
countryAgents.add(ca);
}
@@ -184,21 +186,20 @@ public class CountryAgentManager {
}
- public void serializeTimberAndCarbonProdForAll() {
- Map<CompositeCountry, Double[]> timberAndCarbonProdMap = new HashMap<CompositeCountry, Double[]>();
+ public void serializeWoodUsageForAll() {
+ Map<CompositeCountry, Map<WoodType, WoodUsageData>> woodUsageMap = new HashMap<CompositeCountry, Map<WoodType, WoodUsageData>>();
for (CountryAgent ca : gamsCountryAgents) {
- Double[] timberAndCarbonProd = {ca.getWoodUsageData().getNetImport(), ca.getNetCarbonFlux()};
- timberAndCarbonProdMap.put(ca.country, timberAndCarbonProd);
+ woodUsageMap.put(ca.country, ca.getWoodUsageData());
}
try {
- LogWriter.println("Starting serializing timber and carbon production to " + ModelConfig.SERIALIZED_WOOD_AND_CARBON_PROD_FILE);
- FileOutputStream fileOut = new FileOutputStream(ModelConfig.SERIALIZED_WOOD_AND_CARBON_PROD_FILE);
+ LogWriter.println("Starting serializing WoodUsage to " + ModelConfig.SERIALIZED_WOOD_USAGE_FILE);
+ FileOutputStream fileOut = new FileOutputStream(ModelConfig.SERIALIZED_WOOD_USAGE_FILE);
ObjectOutputStream out = new ObjectOutputStream(fileOut);
- out.writeObject(timberAndCarbonProdMap);
+ out.writeObject(woodUsageMap);
out.close();
fileOut.close();
- LogWriter.println("Serialized timber and carbon production is saved");
+ LogWriter.println("Serialized WoodUsage is saved");
} catch (IOException i) {
i.printStackTrace();
}
diff --git a/src/ac/ed/lurg/country/ForestManager.java b/src/ac/ed/lurg/country/ForestManager.java
deleted file mode 100644
index 834e8063897e49177b2d2f0db2b64806b511b038..0000000000000000000000000000000000000000
--- a/src/ac/ed/lurg/country/ForestManager.java
+++ /dev/null
@@ -1,140 +0,0 @@
-package ac.ed.lurg.country;
-
-import java.util.ArrayList;
-import java.util.HashMap;
-import java.util.Map;
-import java.util.Map.Entry;
-import ac.ed.lurg.ModelConfig;
-import ac.ed.lurg.Timestep;
-import ac.ed.lurg.landuse.LandUseItem;
-import ac.ed.lurg.landuse.WoodYieldRasterSet;
-import ac.ed.lurg.types.LandCoverType;
-import ac.ed.lurg.utils.DoubleMap;
-import ac.ed.lurg.utils.LazyTreeMap;
-import ac.ed.lurg.utils.LogWriter;
-import ac.sac.raster.IntegerRasterItem;
-import ac.sac.raster.RasterKey;
-import ac.sac.raster.RasterSet;
-import ac.ed.lurg.landuse.WoodYieldItem;
-
-public class ForestManager {
- private Map<Integer, Double> timberForestYield;
- private Map<Integer, Double> timberForestArea;
-
- ForestManager() {
- timberForestYield = new HashMap<Integer, Double>();
- timberForestArea = new HashMap<Integer, Double>();
- }
-
- public void initialiseForests(RasterSet<IntegerRasterItem> yieldClusters, RasterSet<LandUseItem> initialLandUse, WoodYieldRasterSet initWoodYieldData) {
- LazyTreeMap<Integer, Double> managedForestArea = new LazyTreeMap<Integer, Double>() {
- private static final long serialVersionUID = 1L;
- protected Double createValue() {return 0.0;}
- };
- LazyTreeMap<Integer, Double> timberStock = new LazyTreeMap<Integer, Double>() {
- private static final long serialVersionUID = 1L;
- protected Double createValue() {return 0.0;}
- };
-
- for (Entry<RasterKey, LandUseItem> entry : initialLandUse.entrySet()) {
- RasterKey key = entry.getKey();
- LandUseItem luItem = entry.getValue();
-
- if (luItem == null) {
- LogWriter.printlnWarning("CountryAgent.initialiseForestManager: null luItem");
- continue;
- }
-
- if (yieldClusters.get(key) != null) {
-
- int clusterId = yieldClusters.get(key).getInt();
- double timberForestAreaSoFar = managedForestArea.lazyGet(clusterId);
- double timberStockSoFar = timberStock.lazyGet(clusterId);
-
- WoodYieldItem wYItem = initWoodYieldData.get(key);
- // TODO Assuming forest was also forest before
- double timberYieldThisTime = (wYItem != null) ? initWoodYieldData.get(key).getWoodYieldRota(LandCoverType.TIMBER_FOREST, LandCoverType.TIMBER_FOREST) : 0.0;
- double timberForestAreaThisTime = luItem.getUnprotectedLandCoverArea(LandCoverType.TIMBER_FOREST);
- managedForestArea.put(clusterId, timberForestAreaSoFar + timberForestAreaThisTime);
- timberStock.put(clusterId, timberStockSoFar + timberForestAreaThisTime * timberYieldThisTime);
-
- } else {
- LogWriter.printlnWarning("CountryAgent.initialiseForestManager: key not mapped to cluster");
- continue;
- }
- }
-
- // Initial timber forest area and yield
- for (Entry<Integer, Double> entry : managedForestArea.entrySet()) {
- int clusterId = entry.getKey();
- double area = entry.getValue();
- if (area > 0) {
- double timberYield = timberStock.get(clusterId) / area;
- timberForestYield.put(clusterId, timberYield);
- timberForestArea.put(clusterId, area);
- }
- }
- }
-
- public DoubleMap<Integer, LandCoverType, Double> getTimberYield() {
- DoubleMap<Integer, LandCoverType, Double> yield = new DoubleMap<Integer, LandCoverType, Double>();
- // Weighted average
- for (Entry<Integer, Double> entry : timberForestYield.entrySet()) {
- yield.put(entry.getKey(), LandCoverType.TIMBER_FOREST, entry.getValue());
- }
-
- return yield;
- }
-
- // New forest yield gradually tends towards base level (based on yield if area was also previously forest)
- private double updateTimberForestYield(double area, double currentYield, double baseYield) {
- double affectedArea = area / ModelConfig.FOREST_ROTATION_PERIOD;
- double unaffectedArea = area - affectedArea;
- return (currentYield * unaffectedArea + baseYield * affectedArea) / area; // weighted average
- }
-
- public void updateForestStatus(ArrayList<LandCoverChangeItem> gamsLandCoverChanges, Map<Integer, Double> baseTimberYield) {
-
- if (!ModelConfig.IS_CALIBRATION_RUN) {
- // Update timber forest yield - the yield effect of previous land use gradually decreases
- for (Integer locId : timberForestYield.keySet()) {
- double area = timberForestArea.get(locId);
- double currentYield = timberForestYield.get(locId);
- double baseYield = baseTimberYield.get(locId);
- double updatedYield = updateTimberForestYield(area, currentYield, baseYield);
- timberForestYield.put(locId, updatedYield);
- }
- }
-
- // Removed timber forest
- for (LandCoverChangeItem item : gamsLandCoverChanges) {
- if (!item.getToLandCover().equals(LandCoverType.TIMBER_FOREST.getName()) && item.getFromLandCover().equals(LandCoverType.TIMBER_FOREST.getName())) {
- int locId = item.getLocation();
- double currentArea = timberForestArea.get(locId);
- double removedArea = item.getArea();
- timberForestArea.put(locId, currentArea - removedArea);
- }
- }
-
- // New timber forest
- for (LandCoverChangeItem item : gamsLandCoverChanges) {
- if (item.getToLandCover().equals(LandCoverType.TIMBER_FOREST.getName()) && !item.getFromLandCover().equals(LandCoverType.TIMBER_FOREST.getName())) {
- int locId = item.getLocation();
- double currentYield = 0;
- double currentArea = 0;
- if (timberForestArea.containsKey(locId)) {
- currentYield = timberForestYield.get(locId);
- currentArea = timberForestArea.get(locId);
- }
- double newArea = item.getArea();
- double newYield = item.getTimberYield();
- double updatedYield = (currentYield * currentArea + newYield * newArea) / (currentArea + newArea); // area weighted average
- timberForestYield.put(locId, updatedYield);
- timberForestArea.put(locId, currentArea + newArea);
- }
- }
-
- // Safety check
- timberForestArea.forEach((k, v) -> {if (v < 0) LogWriter.printlnError("ForestManager.updateForestHistory(): timberForestArea < 0");});
- }
-}
diff --git a/src/ac/ed/lurg/country/crafty/CraftyCountryAgent.java b/src/ac/ed/lurg/country/crafty/CraftyCountryAgent.java
index 2a14d09916880a40110f963b8048f3cdc944ad0b..201cc94b46ba023c40a97aeda832db8520e307e8 100644
--- a/src/ac/ed/lurg/country/crafty/CraftyCountryAgent.java
+++ b/src/ac/ed/lurg/country/crafty/CraftyCountryAgent.java
@@ -12,6 +12,7 @@ import ac.ed.lurg.landuse.CropUsageData;
import ac.ed.lurg.landuse.WoodUsageData;
import ac.ed.lurg.types.CommodityType;
import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.WoodType;
/**
* Country agent that is interface over data transferred to and from CRAFTY.
@@ -57,7 +58,7 @@ public class CraftyCountryAgent extends AbstractCountryAgent {
}
@Override
- public WoodUsageData getWoodUsageData() {
+ public Map<WoodType, WoodUsageData> getWoodUsageData() {
// TODO Auto-generated method stub
return null;
}
diff --git a/src/ac/ed/lurg/country/gams/GamsCountryInput.java b/src/ac/ed/lurg/country/gams/GamsCountryInput.java
index 252b1bbc0d61dbd224a76f1bb43351bc7e788722..7b9f503ca4668eaaf749ff819d7d23a076f9a865 100644
--- a/src/ac/ed/lurg/country/gams/GamsCountryInput.java
+++ b/src/ac/ed/lurg/country/gams/GamsCountryInput.java
@@ -8,6 +8,7 @@ import ac.ed.lurg.country.CompositeCountry;
import ac.ed.lurg.country.CountryPrice;
import ac.ed.lurg.country.TradeConstraint;
import ac.ed.lurg.landuse.CropUsageData;
+import ac.ed.lurg.landuse.WoodUsageData;
import ac.ed.lurg.types.CommodityType;
import ac.ed.lurg.types.CropType;
import ac.ed.lurg.types.WoodType;
@@ -25,13 +26,14 @@ public class GamsCountryInput {
private TradeConstraint carbonTradeConstraint;
private Map<WoodType, Double> woodDemand;
private CountryPrice woodPrice;
- private TradeConstraint timberTradeConstraint;
+ private Map<WoodType, TradeConstraint> woodTradeConstraints;
+ private Map<WoodType, WoodUsageData> previousWoodUsageData;
public GamsCountryInput(CompositeCountry country, Map<CommodityType, Double> projectedDemand, Map<CropType, CountryPrice> countryPrices,
Map<CropType, TradeConstraint> importConstraints, Map<CropType, CropUsageData> previousCropUsageData,
Map<CommodityType, Map<CropType, Double>> minDemandFracts, Map<CropType, Double> subsidyRates,
CountryPrice carbonPrice, TradeConstraint carbonTradeConstraint, Map<WoodType, Double> woodDemand,
- CountryPrice woodPrice, TradeConstraint timberTradeConstraint) {
+ CountryPrice woodPrice, Map<WoodType, TradeConstraint> woodTradeConstraints, Map<WoodType, WoodUsageData> previousWoodUsageData) {
super();
this.country = country;
this.projectedDemand = projectedDemand;
@@ -44,7 +46,8 @@ public class GamsCountryInput {
this.carbonTradeConstraint = carbonTradeConstraint;
this.woodDemand = woodDemand;
this.woodPrice = woodPrice;
- this.timberTradeConstraint = timberTradeConstraint;
+ this.woodTradeConstraints = woodTradeConstraints;
+ this.previousWoodUsageData = previousWoodUsageData;
}
public CompositeCountry getCountry() {
@@ -113,8 +116,12 @@ public class GamsCountryInput {
return woodPrice;
}
- public TradeConstraint getWoodTradeConstraint() {
- return timberTradeConstraint;
+ public Map<WoodType, TradeConstraint> getWoodTradeConstraints() {
+ return woodTradeConstraints;
+ }
+
+ public Map<WoodType, WoodUsageData> getPreviousWoodUsageData() {
+ return previousWoodUsageData;
}
public void setTradeConstraints(Map<CropType, TradeConstraint> tradeConstraints) {
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
index 211385243c1735e39aa1ac2644d46fd5305ba082..6ef2102c49bf8bd48dff4aaad79f2d844064d6ad 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
@@ -341,8 +341,18 @@ public class GamsLocationOptimiser {
addScalar(inDB, "carbonPrice", countryInput.getCarbonPrice().getExportPrice(), 5);
addScalar(inDB, "woodPrice", countryInput.getWoodPrice().getExportPrice(), 5);
- addScalar(inDB, "woodMaxNetImport", countryInput.getWoodTradeConstraint().getMaxConstraint(), 5);
- addScalar(inDB, "woodMinNetImport", countryInput.getWoodTradeConstraint().getMinConstraint(), 5);
+
+ // Not simulating production of different wood products so need to aggregate
+ double woodMinNetImport = 0;
+ double woodMaxNetImport = 0;
+ Map<WoodType, TradeConstraint> woodTradeConstraints = countryInput.getWoodTradeConstraints();
+ for (TradeConstraint tc: woodTradeConstraints.values()) {
+ woodMinNetImport += tc.getMinConstraint();
+ woodMaxNetImport += tc.getMaxConstraint();
+ }
+
+ addScalar(inDB, "woodMaxNetImport", woodMaxNetImport, 5);
+ addScalar(inDB, "woodMinNetImport", woodMinNetImport, 5);
// Wood demand
Map<WoodType, Double> woodDemandMap = countryInput.getWoodDemand();
@@ -424,7 +434,7 @@ public class GamsLocationOptimiser {
conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.CROPLAND, ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.CROP_INCREASE_COST);
conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.PASTURE, ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.PASTURE_INCREASE_COST);
conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.NATURAL, ModelConfig.MANAGED_FOREST_DECREASE_COST);
- conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.TIMBER_FOREST, ModelConfig.FOREST_MANAGEMENT_COST);
+ conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.TIMBER_FOREST, 0.0);
conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.CARBON_FOREST, ModelConfig.FOREST_MANAGEMENT_COST);
conversionCosts.put(LandCoverType.CARBON_FOREST, LandCoverType.CROPLAND, ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.CROP_INCREASE_COST);
@@ -448,6 +458,8 @@ public class GamsLocationOptimiser {
addScalar(inDB, "forestRotationPeriod", ModelConfig.FOREST_ROTATION_PERIOD, 3);
double forestChangeRate = (ModelConfig.IS_CALIBRATION_RUN) ? 0.0 : ModelConfig.FOREST_MAX_CHANGE;
addScalar(inDB, "managedForestAreaMaxChangeRate", forestChangeRate, 3);
+ addScalar(inDB, "forestEstablishmentCost", ModelConfig.FOREST_ESTABLISHMENT_COST, 3);
+ addScalar(inDB, "naturalAreaValue", ModelConfig.NATURAL_AREA_VALUE, 3);
}
private void addScalar(GAMSDatabase gamsDb, String recordName, double val, int places) {
@@ -612,11 +624,27 @@ public class GamsLocationOptimiser {
double netCarbonFlux = outDB.getParameter("netCarbonFlux").getFirstRecord().getValue();
// Timber harvest
- double totalTimberProduction = outDB.getParameter("totalWoodHarvest").getFirstRecord().getValue();
+ double totalWoodHarvest = outDB.getParameter("totalWoodHarvest").getFirstRecord().getValue();
double netWoodImport = outDB.getParameter("netWoodImport").getFirstRecord().getValue();
- WoodUsageData woodUsageData = new WoodUsageData(totalTimberProduction, netWoodImport);
+
+ // Need to disaggregate wood harvest and import by type. This is just for future proofing. Eventually, could have separate harvest for each wood type.
+ Map<WoodType, WoodUsageData> newWoodUsageMap = new HashMap<WoodType, WoodUsageData>();
+ Map<WoodType, WoodUsageData> previousWoodUsageMap = inputData.getCountryInput().getPreviousWoodUsageData();
+ double previousNetImportSum = previousWoodUsageMap.values().stream().mapToDouble(o -> o.getNetImport()).sum();
+ double netWoodImportChange = netWoodImport - previousNetImportSum;
+ double previousWoodHarvestSum = previousWoodUsageMap.values().stream().mapToDouble(o -> o.getHarvest()).sum();
+ int numWoodTypes = WoodType.values().length;
+ for (WoodType wt : WoodType.values()) {
+ double previousNetImport = previousWoodUsageMap.get(wt).getNetImport();
+ double newNetImport = previousNetImport + netWoodImportChange / numWoodTypes;
+ double previousHarvest = previousWoodUsageMap.get(wt).getHarvest();
+ double newHarvest = totalWoodHarvest * (previousHarvest / previousWoodHarvestSum);
+ newWoodUsageMap.put(wt, new WoodUsageData(newHarvest, newNetImport));
+ }
+
+ WoodUsageData woodUsageData = new WoodUsageData(totalWoodHarvest, netWoodImport);
- GamsLocationOutput results = new GamsLocationOutput(modelStatus, landUses, cropUsageData, landCoverChanges, gamsLandCoverChanges, netCarbonFlux, woodUsageData);
+ GamsLocationOutput results = new GamsLocationOutput(modelStatus, landUses, cropUsageData, landCoverChanges, gamsLandCoverChanges, netCarbonFlux, newWoodUsageMap);
return results;
}
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOutput.java b/src/ac/ed/lurg/country/gams/GamsLocationOutput.java
index 0876507f0325a3868d6cede31eaae06fb62d3704..07637c1b33bb45de7991390efc20d87c072c03f3 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationOutput.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationOutput.java
@@ -11,6 +11,7 @@ import ac.ed.lurg.landuse.WoodUsageData;
import ac.ed.lurg.types.CropType;
import ac.ed.lurg.utils.TripleMap;
import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.types.WoodType;
import ac.ed.lurg.country.LandCoverChangeItem;
public class GamsLocationOutput {
@@ -21,14 +22,14 @@ public class GamsLocationOutput {
private TripleMap<Integer, LandCoverType, LandCoverType, Double> landCoverChanges;
ArrayList<LandCoverChangeItem> gamsLandCoverChanges;
private double netCarbonFlux;
- private WoodUsageData woodUsageData;
+ private Map<WoodType, WoodUsageData> woodUsageData;
public GamsLocationOutput(ModelStat status,
Map<Integer, LandUseItem> landUses,
Map<CropType, CropUsageData> cropUsageData,
TripleMap<Integer, LandCoverType, LandCoverType, Double> landCoverChange,
ArrayList<LandCoverChangeItem> gamsLandCoverChanges,
- double netCarbonFlux, WoodUsageData woodUsageData) {
+ double netCarbonFlux, Map<WoodType, WoodUsageData> woodUsageData) {
super();
this.status = status;
this.landUses = landUses;
@@ -62,7 +63,7 @@ public class GamsLocationOutput {
return netCarbonFlux;
}
- public WoodUsageData getWoodUsageData() {
+ public Map<WoodType, WoodUsageData> getWoodUsageData() {
return woodUsageData;
}
}
diff --git a/src/ac/ed/lurg/country/gams/GamsRasterOutput.java b/src/ac/ed/lurg/country/gams/GamsRasterOutput.java
index cade1516dc2e41e9fab56c3634f3b88145617bbf..0add17f98234b42295ce9211b31726c595ef873c 100644
--- a/src/ac/ed/lurg/country/gams/GamsRasterOutput.java
+++ b/src/ac/ed/lurg/country/gams/GamsRasterOutput.java
@@ -10,6 +10,7 @@ import ac.ed.lurg.landuse.CropUsageData;
import ac.ed.lurg.landuse.LandUseItem;
import ac.ed.lurg.landuse.WoodUsageData;
import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.WoodType;
import ac.sac.raster.RasterSet;
public class GamsRasterOutput {
@@ -20,28 +21,28 @@ public class GamsRasterOutput {
private ArrayList<LandCoverChangeItem> gamsLandCoverChanges;
private Map<Integer, Double> baseTimberYield;
private double netCarbonFlux;
- private WoodUsageData woodUsageData;
+ private Map<WoodType, WoodUsageData> woodUsageData;
- public GamsRasterOutput(RasterSet<LandUseItem> landUses, Map<CropType, CropUsageData> cropUsageData) {
+ public GamsRasterOutput(RasterSet<LandUseItem> landUses, Map<CropType, CropUsageData> cropUsageData, Map<WoodType, WoodUsageData> woodUsageData) {
super();
this.landUses = landUses;
this.cropUsageData = cropUsageData;
+ this.woodUsageData = woodUsageData;
}
-
+/*
public GamsRasterOutput(RasterSet<LandUseItem> landUses, Map<CropType, CropUsageData> cropUsageData, Double[] woodAndCarbonProd) {
this(landUses, cropUsageData);
this.woodUsageData = new WoodUsageData(woodAndCarbonProd[0], 0.0);
this.netCarbonFlux = woodAndCarbonProd[1];
}
-
+*/
public GamsRasterOutput(ModelStat status, RasterSet<LandUseItem> intensityRaster, Map<CropType, CropUsageData> cropUsageData,
- ArrayList<LandCoverChangeItem> gamsLandCoverChanges, Map<Integer, Double> baseTimberYield, double netCarbonFlux, WoodUsageData woodUsageData) {
- this(intensityRaster, cropUsageData);
+ ArrayList<LandCoverChangeItem> gamsLandCoverChanges, Map<Integer, Double> baseTimberYield, double netCarbonFlux, Map<WoodType, WoodUsageData> woodUsageData) {
+ this(intensityRaster, cropUsageData, woodUsageData);
this.status = status;
this.gamsLandCoverChanges = gamsLandCoverChanges;
this.baseTimberYield = baseTimberYield;
this.netCarbonFlux = netCarbonFlux;
- this.woodUsageData = woodUsageData;
}
public ModelStat getStatus() {
@@ -68,7 +69,7 @@ public class GamsRasterOutput {
return netCarbonFlux;
}
- public WoodUsageData getWoodUsageData() {
+ public Map<WoodType, WoodUsageData> getWoodUsageData() {
return woodUsageData;
}
diff --git a/src/ac/ed/lurg/demand/AbstractDemandManager.java b/src/ac/ed/lurg/demand/AbstractDemandManager.java
index 99f8206a5272024eff35083c73e944f16c6f7651..5367ff3086305bd566fe9d6caf768e05f207b8d6 100644
--- a/src/ac/ed/lurg/demand/AbstractDemandManager.java
+++ b/src/ac/ed/lurg/demand/AbstractDemandManager.java
@@ -151,6 +151,10 @@ public abstract class AbstractDemandManager {
for (SingleCountry c : compositeCountryManager.getAllForCompositeCountry(cc)) {
singleDemandMap = getWoodDemand(c, year);
+ if (singleDemandMap.isEmpty()) {
+ LogWriter.printlnWarning("Wood demand not found for: " + c.getCountryName());
+ continue;
+ }
for (WoodType w : WoodType.values()) {
double totalDemand = (compositeDemandMap.containsKey(w)) ? compositeDemandMap.get(w) : 0.0;
double demand = singleDemandMap.get(w);
diff --git a/src/ac/ed/lurg/demand/WoodDemandManager.java b/src/ac/ed/lurg/demand/WoodDemandManager.java
index b95093a3a91b6191920da7fb47e5203067732031..4a2caae2245dc2097adf09921ca1c4715ae8e9a5 100644
--- a/src/ac/ed/lurg/demand/WoodDemandManager.java
+++ b/src/ac/ed/lurg/demand/WoodDemandManager.java
@@ -16,7 +16,6 @@ public class WoodDemandManager {
private static final int COUNTY_COL = 0;
private static final int WOOD_TYPE_COL = 1;
private static final int DEMAND_COL = 2;
- private static final double WOOD_BIOMASS_CONVERSION_FACTOR = 0.2688; // m3 to tC-eq p.16 [https://doi.org/10.5194/gmd-13-5425-2020]
private Map<SingleCountry, Map<WoodType, Double>> baseWoodDemand = new HashMap<SingleCountry, Map<WoodType, Double>>();
@@ -55,7 +54,7 @@ public class WoodDemandManager {
baseWoodDemand.put(country, woodTypeMap);
}
- double demandBiomassEq = demand / 1e6 * WOOD_BIOMASS_CONVERSION_FACTOR; // converting from m3 to MtC-eq
+ double demandBiomassEq = demand * ModelConfig.WOOD_BIOMASS_CONVERSION_FACTOR; // converting from m3 to MtC-eq
woodTypeMap.put(woodType, demandBiomassEq);
}
}
diff --git a/src/ac/ed/lurg/forestry/CarbonForestGrowthDataReader.java b/src/ac/ed/lurg/forestry/CarbonForestGrowthDataReader.java
new file mode 100644
index 0000000000000000000000000000000000000000..8b65ec9df73137c13146c4edaa04912c5653547a
--- /dev/null
+++ b/src/ac/ed/lurg/forestry/CarbonForestGrowthDataReader.java
@@ -0,0 +1,43 @@
+package ac.ed.lurg.forestry;
+
+import java.util.Map;
+
+import ac.ed.lurg.types.LandCoverType;
+import ac.sac.raster.AbstractTabularRasterReader;
+import ac.sac.raster.RasterKey;
+import ac.sac.raster.RasterSet;
+
+public class CarbonForestGrowthDataReader extends AbstractTabularRasterReader<ForestGrowthItem>{
+ private static final int MIN_COLS = 22;
+ private static final double conversionFactor = 10.0; // convert kgC/m2 to tC/ha
+
+ public CarbonForestGrowthDataReader(RasterSet<ForestGrowthItem> forestGrowthData) {
+ super("[ |\t]+", MIN_COLS, forestGrowthData);
+ }
+
+ protected void setData(RasterKey key, ForestGrowthItem item, Map<String, Double> rowValues) {
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 0, 0.0);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 5, getValueForCol(rowValues, "5") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 10, getValueForCol(rowValues, "10") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 15, getValueForCol(rowValues, "15") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 20, getValueForCol(rowValues, "20") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 25, getValueForCol(rowValues, "25") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 30, getValueForCol(rowValues, "30") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 35, getValueForCol(rowValues, "35") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 40, getValueForCol(rowValues, "40") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 45, getValueForCol(rowValues, "45") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 50, getValueForCol(rowValues, "50") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 55, getValueForCol(rowValues, "55") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 60, getValueForCol(rowValues, "60") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 65, getValueForCol(rowValues, "65") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 70, getValueForCol(rowValues, "70") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 75, getValueForCol(rowValues, "75") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 80, getValueForCol(rowValues, "80") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 85, getValueForCol(rowValues, "85") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 90, getValueForCol(rowValues, "90") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 95, getValueForCol(rowValues, "95") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 100, getValueForCol(rowValues, "100") * conversionFactor);
+
+
+ }
+}
diff --git a/src/ac/ed/lurg/forestry/ForestGrowthItem.java b/src/ac/ed/lurg/forestry/ForestGrowthItem.java
new file mode 100644
index 0000000000000000000000000000000000000000..97cdfbf5fdd4a6dab1b4bb36c83a81afa52996b6
--- /dev/null
+++ b/src/ac/ed/lurg/forestry/ForestGrowthItem.java
@@ -0,0 +1,32 @@
+package ac.ed.lurg.forestry;
+
+import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.utils.DoubleMap;
+import ac.ed.lurg.utils.Interpolator;
+import ac.sac.raster.RasterItem;
+
+public class ForestGrowthItem implements RasterItem {
+ private DoubleMap<LandCoverType, Integer, Double> forestGrowthFunction = new DoubleMap<LandCoverType, Integer, Double>();
+
+ public void setForestGrowthFunction(LandCoverType forestType, int time, double deltaC) {
+ forestGrowthFunction.put(forestType, time, deltaC);
+ }
+
+ public double getGrowth(LandCoverType forestType, int time) {
+ int lowerTime = (time/5) * 5;
+ int upperTime = lowerTime + 5;
+ Double lowerVal = forestGrowthFunction.get(forestType, lowerTime);
+ Double upperVal = forestGrowthFunction.get(forestType, upperTime);
+ double factor = ((double)(time - lowerTime)) / (upperTime - lowerTime);
+ Double d = Interpolator.interpolate(lowerVal, upperVal, factor);
+ return d;
+ }
+
+ public double getCumulGrowth(LandCoverType forestType, int startTime, int endTime) {
+ double cumulGrowth = 0;
+ for (int i=startTime; i<=endTime; i++) {
+ cumulGrowth += getGrowth(forestType, i);
+ }
+ return cumulGrowth;
+ }
+}
diff --git a/src/ac/ed/lurg/forestry/ForestGrowthRasterSet.java b/src/ac/ed/lurg/forestry/ForestGrowthRasterSet.java
new file mode 100644
index 0000000000000000000000000000000000000000..d66c4c5518d8f46cae9bf800925d075c1fa3bf78
--- /dev/null
+++ b/src/ac/ed/lurg/forestry/ForestGrowthRasterSet.java
@@ -0,0 +1,27 @@
+package ac.ed.lurg.forestry;
+
+import java.util.Collection;
+
+import ac.sac.raster.RasterHeaderDetails;
+import ac.sac.raster.RasterKey;
+import ac.sac.raster.RasterSet;
+
+public class ForestGrowthRasterSet extends RasterSet<ForestGrowthItem> {
+
+ private static final long serialVersionUID = -3816941524974512722L;
+
+ public ForestGrowthRasterSet(RasterHeaderDetails header) {
+ super(header);
+ }
+
+ protected ForestGrowthItem createRasterData() {
+ return new ForestGrowthItem();
+ }
+
+ @Override
+ public ForestGrowthRasterSet createSubsetForKeys(Collection<RasterKey> keys) {
+ ForestGrowthRasterSet subsetForestGrowthRaster = new ForestGrowthRasterSet(getHeaderDetails());
+ popSubsetForKeys(subsetForestGrowthRaster, keys);
+ return subsetForestGrowthRaster;
+ }
+}
diff --git a/src/ac/ed/lurg/forestry/ForestManager.java b/src/ac/ed/lurg/forestry/ForestManager.java
new file mode 100644
index 0000000000000000000000000000000000000000..a8ac616c9f6d5b0bc55e30194b617d77e72f17ec
--- /dev/null
+++ b/src/ac/ed/lurg/forestry/ForestManager.java
@@ -0,0 +1,65 @@
+package ac.ed.lurg.forestry;
+
+import java.util.ArrayList;
+import java.util.HashMap;
+import java.util.List;
+import java.util.Map;
+import java.util.Map.Entry;
+import ac.ed.lurg.ModelConfig;
+import ac.ed.lurg.Timestep;
+import ac.ed.lurg.landuse.LandCoverChangeItem;
+import ac.ed.lurg.landuse.LandUseItem;
+import ac.ed.lurg.landuse.WoodYieldItem;
+import ac.ed.lurg.landuse.WoodYieldRasterSet;
+import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.utils.DoubleMap;
+import ac.ed.lurg.utils.LazyTreeMap;
+import ac.ed.lurg.utils.LogWriter;
+import ac.sac.raster.IntegerRasterItem;
+import ac.sac.raster.RasterKey;
+import ac.sac.raster.RasterSet;
+
+public class ForestManager {
+ private Map<LandCoverType, ArrayList<ForestStandItem>> forest;
+ private ForestGrowthItem growthFunction;
+
+ ForestManager(Map<LandCoverType, Double> landCoverAreas, double initialForestedNaturalArea, ForestGrowthItem growthFunction) {
+ forest = new HashMap<LandCoverType, ArrayList<ForestStandItem>>();
+ initialiseForests(landCoverAreas, initialForestedNaturalArea);
+ }
+
+ public void initialiseForests(Map<LandCoverType, Double> landCoverAreas, double initialForestedNaturalArea) {
+
+ for (LandCoverType forestType : LandCoverType.getManagedForestTypes()) {
+ ArrayList<ForestStandItem> forestStands = new ArrayList<ForestStandItem>();
+ double yield = growthFunction.getCumulGrowth(forestType, 0, ModelConfig.FOREST_INIT_AGE);
+ forestStands.add(new ForestStandItem(forestType, landCoverAreas.get(forestType), ModelConfig.FOREST_INIT_AGE, yield));
+ forest.put(forestType, forestStands);
+ }
+
+ ArrayList<ForestStandItem> naturalStands= new ArrayList<ForestStandItem>();
+ double natForestedArea = initialForestedNaturalArea;
+ double natOtherArea = landCoverAreas.get(LandCoverType.NATURAL) - initialForestedNaturalArea;
+
+ double yield = growthFunction.getCumulGrowth(LandCoverType.NATURAL, 0, ModelConfig.FOREST_INIT_AGE);
+ naturalStands.add(new ForestStandItem(LandCoverType.NATURAL, natForestedArea, ModelConfig.FOREST_INIT_AGE, yield));
+
+ naturalStands.add(new ForestStandItem(LandCoverType.NATURAL, natOtherArea, 0, 0.0));
+
+ forest.put(LandCoverType.NATURAL, naturalStands);
+ }
+
+ public void growForests() {
+ for (LandCoverType forestType : LandCoverType.getNaturalTypes()) {
+ ArrayList<ForestStandItem> forestStands = forest.get(forestType);
+ for (ForestStandItem stand: forestStands) {
+ stand.updateAgeAndYield(growthFunction);
+ }
+ }
+ }
+
+ public void setGrowthFunction(ForestGrowthItem growthFunction) {
+ this.growthFunction = growthFunction;
+ }
+
+}
diff --git a/src/ac/ed/lurg/forestry/ForestStandItem.java b/src/ac/ed/lurg/forestry/ForestStandItem.java
new file mode 100644
index 0000000000000000000000000000000000000000..f27b9e7667264d7fc864fbd0ac03cccf1c1808e5
--- /dev/null
+++ b/src/ac/ed/lurg/forestry/ForestStandItem.java
@@ -0,0 +1,51 @@
+package ac.ed.lurg.forestry;
+
+import ac.ed.lurg.ModelConfig;
+import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.utils.LogWriter;
+
+public class ForestStandItem {
+ private LandCoverType forestType;
+ private double area;
+ private int age;
+ private double yield;
+
+ public ForestStandItem(LandCoverType forestType, double area, int age, double yield) {
+ this.forestType = forestType;
+ this.area = area;
+ this.age = age;
+ this.yield = yield;
+ }
+
+ public LandCoverType getForestType() {
+ return forestType;
+ }
+
+ public double getArea() {
+ return area;
+ }
+
+ public double getAge() {
+ return age;
+ }
+
+ public double getYield() {
+ return yield;
+ }
+
+ public void removeArea(double a) {
+ if (a > area) {
+ area -= Math.min(area, a);
+ LogWriter.printlnError("ForestItem: attempting to remove too much area");
+ } else {
+ area -= a;
+ }
+ }
+
+ public void updateAgeAndYield(ForestGrowthItem forestGrowthFunction) {
+ age += ModelConfig.TIMESTEP_SIZE;
+ yield += forestGrowthFunction.getGrowth(forestType, age);
+
+ }
+
+}
diff --git a/src/ac/ed/lurg/forestry/NaturalForestGrowthDataReader.java b/src/ac/ed/lurg/forestry/NaturalForestGrowthDataReader.java
new file mode 100644
index 0000000000000000000000000000000000000000..f2de87ce9e96c0457ef57ed64a3791ed669998bd
--- /dev/null
+++ b/src/ac/ed/lurg/forestry/NaturalForestGrowthDataReader.java
@@ -0,0 +1,43 @@
+package ac.ed.lurg.forestry;
+
+import java.util.Map;
+
+import ac.ed.lurg.types.LandCoverType;
+import ac.sac.raster.AbstractTabularRasterReader;
+import ac.sac.raster.RasterKey;
+import ac.sac.raster.RasterSet;
+
+public class NaturalForestGrowthDataReader extends AbstractTabularRasterReader<ForestGrowthItem>{
+ private static final int MIN_COLS = 22;
+ private static final double conversionFactor = 10.0; // convert kgC/m2 to tC/ha
+
+ public NaturalForestGrowthDataReader(RasterSet<ForestGrowthItem> forestGrowthData) {
+ super("[ |\t]+", MIN_COLS, forestGrowthData);
+ }
+
+ protected void setData(RasterKey key, ForestGrowthItem item, Map<String, Double> rowValues) {
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 0, 0.0);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 5, getValueForCol(rowValues, "5") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 10, getValueForCol(rowValues, "10") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 15, getValueForCol(rowValues, "15") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 20, getValueForCol(rowValues, "20") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 25, getValueForCol(rowValues, "25") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 30, getValueForCol(rowValues, "30") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 35, getValueForCol(rowValues, "35") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 40, getValueForCol(rowValues, "40") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 45, getValueForCol(rowValues, "45") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 50, getValueForCol(rowValues, "50") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 55, getValueForCol(rowValues, "55") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 60, getValueForCol(rowValues, "60") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 65, getValueForCol(rowValues, "65") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 70, getValueForCol(rowValues, "70") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 75, getValueForCol(rowValues, "75") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 80, getValueForCol(rowValues, "80") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 85, getValueForCol(rowValues, "85") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 90, getValueForCol(rowValues, "90") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 95, getValueForCol(rowValues, "95") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.NATURAL, 100, getValueForCol(rowValues, "100") * conversionFactor);
+
+
+ }
+}
diff --git a/src/ac/ed/lurg/forestry/TimberForestGrowthDataReader.java b/src/ac/ed/lurg/forestry/TimberForestGrowthDataReader.java
new file mode 100644
index 0000000000000000000000000000000000000000..f2d8135d474dd3290840ec9c1a92fa55254795ce
--- /dev/null
+++ b/src/ac/ed/lurg/forestry/TimberForestGrowthDataReader.java
@@ -0,0 +1,43 @@
+package ac.ed.lurg.forestry;
+
+import java.util.Map;
+
+import ac.ed.lurg.types.LandCoverType;
+import ac.sac.raster.AbstractTabularRasterReader;
+import ac.sac.raster.RasterKey;
+import ac.sac.raster.RasterSet;
+
+public class TimberForestGrowthDataReader extends AbstractTabularRasterReader<ForestGrowthItem>{
+ private static final int MIN_COLS = 22;
+ private static final double conversionFactor = 10.0; // convert kgC/m2 to tC/ha
+
+ public TimberForestGrowthDataReader(RasterSet<ForestGrowthItem> forestGrowthData) {
+ super("[ |\t]+", MIN_COLS, forestGrowthData);
+ }
+
+ protected void setData(RasterKey key, ForestGrowthItem item, Map<String, Double> rowValues) {
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 0, 0.0);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 5, getValueForCol(rowValues, "5") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 10, getValueForCol(rowValues, "10") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 15, getValueForCol(rowValues, "15") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 20, getValueForCol(rowValues, "20") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 25, getValueForCol(rowValues, "25") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 30, getValueForCol(rowValues, "30") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 35, getValueForCol(rowValues, "35") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 40, getValueForCol(rowValues, "40") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 45, getValueForCol(rowValues, "45") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 50, getValueForCol(rowValues, "50") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 55, getValueForCol(rowValues, "55") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 60, getValueForCol(rowValues, "60") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 65, getValueForCol(rowValues, "65") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 70, getValueForCol(rowValues, "70") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 75, getValueForCol(rowValues, "75") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 80, getValueForCol(rowValues, "80") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 85, getValueForCol(rowValues, "85") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 90, getValueForCol(rowValues, "90") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 95, getValueForCol(rowValues, "95") * conversionFactor);
+ item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 100, getValueForCol(rowValues, "100") * conversionFactor);
+
+
+ }
+}
diff --git a/src/ac/ed/lurg/landuse/CropUsageReader.java b/src/ac/ed/lurg/landuse/CropUsageReader.java
index 195475162e0fe17cf77c44cfde047dd75a2aa9f4..8c243f812b4a83da9ef3f19a62790e873b225a5c 100644
--- a/src/ac/ed/lurg/landuse/CropUsageReader.java
+++ b/src/ac/ed/lurg/landuse/CropUsageReader.java
@@ -76,7 +76,7 @@ public class CropUsageReader {
fitReader.close();
} catch (IOException e) {
- LogWriter.printlnError("Failed in reading commodity demand fits");
+ LogWriter.printlnError("Failed in reading commodity usage data");
LogWriter.print(e);
}
LogWriter.println("Processed " + filename + ", create " + commodityMap.size() + " country commodity maps values");
diff --git a/src/ac/ed/lurg/landuse/LandUseItem.java b/src/ac/ed/lurg/landuse/LandUseItem.java
index ec47344813815d3aef4d50cf908a134d8eede8f1..59372649ad0030610b1c8612d0887ee248104f0a 100644
--- a/src/ac/ed/lurg/landuse/LandUseItem.java
+++ b/src/ac/ed/lurg/landuse/LandUseItem.java
@@ -5,6 +5,7 @@ import java.util.Collection;
import java.util.HashMap;
import java.util.Map;
import ac.ed.lurg.ModelConfig;
+import ac.ed.lurg.forestry.ForestManager;
import ac.ed.lurg.types.CropToDouble;
import ac.ed.lurg.types.CropType;
import ac.ed.lurg.types.LandCoverType;
@@ -23,6 +24,7 @@ public class LandUseItem implements InterpolatingRasterItem<LandUseItem>, Serial
private double unavailableArea; //area unavailable due to altitude etc
private double suitableArea; //Mha
private double initialForestedNaturalArea;
+ private ForestManager forestManager;
public LandUseItem() {}
diff --git a/src/ac/ed/lurg/landuse/WoodUsageData.java b/src/ac/ed/lurg/landuse/WoodUsageData.java
index 2232acd97c23aa74aad946248fa87bbff2a48eb3..2a2d29f86c48b5f3910998b3efcf7fe2c3a98276 100644
--- a/src/ac/ed/lurg/landuse/WoodUsageData.java
+++ b/src/ac/ed/lurg/landuse/WoodUsageData.java
@@ -1,6 +1,11 @@
package ac.ed.lurg.landuse;
-public class WoodUsageData {
+import java.io.Serializable;
+
+public class WoodUsageData implements Serializable {
+
+ private static final long serialVersionUID = -3329080279189782763L;
+
private double harvest;
private double netImport;
diff --git a/src/ac/ed/lurg/landuse/WoodUsageReader.java b/src/ac/ed/lurg/landuse/WoodUsageReader.java
new file mode 100644
index 0000000000000000000000000000000000000000..3652e3c6b240c778407e4d8ed0c03f66f2dc1771
--- /dev/null
+++ b/src/ac/ed/lurg/landuse/WoodUsageReader.java
@@ -0,0 +1,94 @@
+package ac.ed.lurg.landuse;
+
+import java.io.BufferedReader;
+import java.io.FileReader;
+import java.io.IOException;
+import java.util.HashMap;
+import java.util.Map;
+
+import ac.ed.lurg.country.CompositeCountryManager;
+import ac.ed.lurg.country.CountryManager;
+import ac.ed.lurg.country.SingleCountry;
+import ac.ed.lurg.ModelConfig;
+import ac.ed.lurg.country.CompositeCountry;
+import ac.ed.lurg.types.WoodType;
+import ac.ed.lurg.utils.LazyHashMap;
+import ac.ed.lurg.utils.LogWriter;
+
+public class WoodUsageReader {
+ private static final int COUNTRY_COL = 0;
+ private static final int WOOD_TYPE_COL = 1;
+ private static final int HARVEST_COL = 2;
+ private static final int NET_IMPORT_COL = 3;
+
+ private CompositeCountryManager compositeCountryManager;
+
+ public WoodUsageReader(CompositeCountryManager compositeCountryManager) {
+ this.compositeCountryManager = compositeCountryManager;
+ }
+
+ public Map<CompositeCountry, Map<WoodType, WoodUsageData>> getWoodUsageData() {
+
+ LazyHashMap<CompositeCountry, Map<WoodType, WoodUsageData>> usageMap = new LazyHashMap<CompositeCountry, Map<WoodType, WoodUsageData>>() {
+ protected Map<WoodType, WoodUsageData> createValue() {
+ Map<WoodType, WoodUsageData> countryData = new HashMap<WoodType, WoodUsageData>();
+ return countryData;
+ }
+ };
+
+ String filename = ModelConfig.WOOD_NET_IMPORTS_FILE;
+ try {
+ BufferedReader reader = new BufferedReader(new FileReader(filename));
+ String line, countryName, woodTypeName;
+ double harvest, netImport;
+ reader.readLine(); // read header
+
+ while ((line = reader.readLine()) != null) {
+ String[] tokens = line.split(",");
+
+ if (tokens.length < 4)
+ LogWriter.printlnError("Too few columns in file: " + filename + ", line: " + line);
+
+ countryName = tokens[COUNTRY_COL];
+ woodTypeName = tokens[WOOD_TYPE_COL];
+ harvest = Double.parseDouble(tokens[HARVEST_COL]) * ModelConfig.WOOD_BIOMASS_CONVERSION_FACTOR; // m3 to MtC-eq
+ netImport = Double.parseDouble(tokens[NET_IMPORT_COL]) * ModelConfig.WOOD_BIOMASS_CONVERSION_FACTOR; // m3 to MtC-eq
+
+ SingleCountry country = CountryManager.getForName(countryName);
+ WoodType woodType = WoodType.getForName(woodTypeName);
+
+ if (country == null) {
+ LogWriter.printlnWarning("WoodUsageReader can't find single country: " + countryName);
+ continue;
+ }
+
+ CompositeCountry cc = compositeCountryManager.getForSingleCountry(country);
+
+ if (cc == null) {
+ LogWriter.printlnWarning("WoodUsageReader can't find composite country for: " + country.getCountryName());
+ continue;
+ }
+
+ Map<WoodType, WoodUsageData> countryData = usageMap.lazyGet(cc);
+ WoodUsageData oldData = countryData.get(woodType);
+
+ // aggregate if multiple countries for cc
+ if (oldData != null) {
+ netImport += oldData.getNetImport();
+ harvest += oldData.getHarvest();
+ }
+
+ WoodUsageData newData = new WoodUsageData(harvest, netImport);
+ countryData.put(woodType, newData);
+ }
+ reader.close();
+
+ } catch (IOException e) {
+ LogWriter.printlnError("Failed in reading wood usage data");
+ LogWriter.print(e);
+ }
+ LogWriter.println("Processed " + filename + ", create " + usageMap.size() + " country wood usage map values");
+
+ return usageMap;
+ }
+}
diff --git a/src/ac/ed/lurg/types/LandCoverType.java b/src/ac/ed/lurg/types/LandCoverType.java
index 297f1dc9b4e33a33f12099fe649cdf5283f74a6d..fc9346e8212e6f212e46d0e8dfdac396c7e7df1e 100644
--- a/src/ac/ed/lurg/types/LandCoverType.java
+++ b/src/ac/ed/lurg/types/LandCoverType.java
@@ -100,5 +100,17 @@ public enum LandCoverType {
return managedForestTypes;
}
+
+ public static Collection<LandCoverType> getNaturalTypes() {
+
+ Collection<LandCoverType> naturalTypes = new HashSet<LandCoverType>();
+
+ for (LandCoverType c : values())
+ if (c.isNatural)
+ naturalTypes.add(c);
+
+ return naturalTypes;
+
+ }
}