lemon-project-template-glpk

comparison deps/glpk/examples/dbf/ForestMgt_Model_I_GIS_dbf.mod @ 10:5545663ca997

find changesets by author, revision, files, or words in the commit message
Configure GLPK build
author Alpar Juttner <alpar@cs.elte.hu>
date Sun, 06 Nov 2011 21:42:23 +0100 (2011-11-06)
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:ab386d7cc2c8
1 # Model I Forest Estate Modelling using GLPK/MathProg
2 # Reading and writing dbf files
3
4 # by Noli Sicad --- nsicad@gmail.com
5 # 18 December 2009
6
7 # Forest Management 4th Edition
8 # by Lawrence Davis, K. Norman Johnson, Pete Bettinger, Theodore Howard
9 # Chapter 11 - Daniel Pickett
10 # http://warnell.forestry.uga.edu/Warnell/Bettinger/mgtbook/index.htm
11
12 # Model I Formulation
13
14 /* Note: This is not the full LP model mentioned in the book.
15 Some of the constraints are deliberately omitted in this model for the purpose of clarity.
16
17 The features of MathProg in this example are:
18 * reading and writing dbf from regular dbf files,
19 * reading dbf file (database of shapefile (stands.shp)) (e.g. area parameter),
20 * using the area data in the constraints and
21 * writing dbf file from result of LP model.
22
23 Model I - Harvest Scheduling formulation for Sustainable Forest Management (SFM)
24
25 Features are:
26 * Net Present Value for the objective function (Revenue - Cost)
27 * Harvest Constraints by period - Sustainable Yield per Period
28 * Even-Flow Constraint / Volume - Harvest Flow Constraint - Alpha (1-Apha)
29 * Even-Flow Constraint / Volume - Harvest Flow Constraint - Beta (1 +Beta)
30 * Forest / Land Constraint -- Total Area of the forest
31 * Forest Stand Constraint -- Individual stands
32
33 What is next? -- Forest Mgt Carbon Accounting for Climate Change
34
35 Note: The model file that the data containing in
36 the dbf files is public domain material (so it is compatible with the
37 GNU GPL) and data can be found in
38 http://warnell.forestry.uga.edu/Warnell/Bettinger/mgtbook/index.htm
39
40 # Noli Sicad --- nsicad@gmail.com
41
42 */
43
44 set G_STAND_TYPE; # A, B, C
45
46 set I_CULTURAL_PRES;
47 set J_MGT_YEAR;
48
49 param K_PERIOD;
50 param Forest_Cost{G_STAND_TYPE,I_CULTURAL_PRES, J_MGT_YEAR}; # cost
51
52 param Yield_Table_Vol{G_STAND_TYPE, I_CULTURAL_PRES, J_MGT_YEAR, 1..K_PERIOD} >= 0;
53
54
55 param Alpha >= 0;
56 param Beta >= 0;
57
58 param TCost_Table{G_STAND_TYPE, I_CULTURAL_PRES, J_MGT_YEAR, 1..K_PERIOD} >= 0;
59
60 param NetRev_Table{G_STAND_TYPE, I_CULTURAL_PRES, J_MGT_YEAR, 1..K_PERIOD};
61
62
63 var XForestLand{g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} >= 0;
64
65
66 #reading dbf tables
67 table tab IN "xBASE" "standtype.dbf": G_STAND_TYPE <- [STAND];
68 display G_STAND_TYPE;
69
70
71 table tab2 IN "xBASE" "cultural_pres.dbf": I_CULTURAL_PRES <- [CUL_PRES];
72 display I_CULTURAL_PRES;
73
74 table tab3 IN "xBASE" "mgt_year.dbf": J_MGT_YEAR <- [MGT_YEAR];
75 display J_MGT_YEAR;
76
77 /*
78 param Forest_Cost{G_STAND_TYPE,I_CULTURAL_PRES, J_MGT_YEAR} default 0; # cost
79 */
80
81 set S1, dimen 3;
82 table tab4 IN "xBASE" "Forest_Cost.dbf": S1 <- [STAND, CUL_PRES, MGT_YEAR],Forest_Cost ~FCOST;
83 display Forest_Cost;
84
85 set S2, dimen 4;
86 table tab5 IN "xBASE" "Yield_Table_Vol.dbf": S2 <- [STAND, CUL_PRES, MGT_YEAR, PERIOD],Yield_Table_Vol ~YIELD;
87 display Yield_Table_Vol;
88
89 set S3, dimen 4;
90 table tab5 IN "xBASE" "TCost_Table.dbf": S3 <- [STAND, CUL_PRES, MGT_YEAR, PERIOD],TCost_Table ~TCOST;
91 display TCost_Table;
92
93
94 set S4, dimen 4;
95 table tab5 IN "xBASE" "NetRev_Table.dbf": S4 <- [STAND, CUL_PRES, MGT_YEAR, PERIOD],NetRev_Table ~NETREV;
96 display NetRev_Table;
97
98
99 param MGT;
100
101 param Area_Stand_Indi{g in G_STAND_TYPE, m in 1..MGT} default 0;
102
103 set ST, dimen 2;
104 table tab5 IN "xBASE" "stands.dbf": ST <- [VEG_TYPE, MGT], Area_Stand_Indi ~ACRES;
105 display Area_Stand_Indi;
106
107 param Area_Stand_Type{g in G_STAND_TYPE}:= sum {m in 1..MGT } Area_Stand_Indi[g,m];
108 display Area_Stand_Type;
109
110
111 param Total_Area := sum {g in G_STAND_TYPE, m in 1..MGT } Area_Stand_Indi[g,m];
112 display Total_Area;
113
114 param Harvest_Min_Vol_Period;
115
116
117 var NetPresentValue;
118
119 # Objective function
120 maximize Net_Present_Value: NetPresentValue;
121
122 subject to NPV:
123 NetPresentValue = sum {g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} Forest_Cost[g,i,j] * XForestLand[g,i,j];
124
125 # Harvest Constraint by Period
126 subject to Harvest_Period_H {k in 1..K_PERIOD}:
127 sum {g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} Yield_Table_Vol[g,i,j,k] * XForestLand[g,i,j] >= Harvest_Min_Vol_Period;
128
129
130 #Even-Flow Constraint / Volume - Harvest Flow Constraint - Alpha
131 subject to Even_Flow_Constaints_Alpha {k in 6..K_PERIOD-1}:
132 (1 - Alpha) * sum {g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} Yield_Table_Vol[g,i,j,k] * XForestLand[g,i,j] -
133 sum {g in G_STAND_TYPE,i in I_CULTURAL_PRES, j in J_MGT_YEAR} Yield_Table_Vol[g,i,j,k+1] * XForestLand[g,i,j] <= 0;
134
135 # Even-Flow Constraint / Volume - Harvest Flow Constraint - Beta
136 subject to Even_Flow_Constaints_Beta {k in 6..K_PERIOD-1}:
137 (1 + Beta) * sum {g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} Yield_Table_Vol[g,i,j,k] * XForestLand[g,i,j] -
138 sum {g in G_STAND_TYPE,i in I_CULTURAL_PRES, j in J_MGT_YEAR} Yield_Table_Vol[g,i,j,k+1] * XForestLand[g,i,j] >= 0;
139
140 # Forest / Land Constraints
141 subject to Total_Area_Constraint:
142 sum {g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} XForestLand[g,i,j] <= Total_Area;
143 display Total_Area;
144
145 # Forest / Land Constraints for A B C
146 subject to Area {g in G_STAND_TYPE}:
147 sum {i in I_CULTURAL_PRES,j in J_MGT_YEAR} XForestLand[g,i,j] = Area_Stand_Type[g];
148
149
150
151 solve;
152 #RESULT SECTION
153 printf '#################################\n';
154 printf 'Forest Management Model I - Noli Sicad\n';
155 printf '\n';
156 printf 'Net Present Value = %.2f\n', NetPresentValue;
157 printf '\n';
158
159 printf '\n';
160 printf 'Variables\n';
161 printf 'Stand_Type Age_Class Mgt_Presc Sign Value \n';
162 printf{g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR}:'%5s %10s %11s = %10.2f\n', g,i,j, XForestLand[g,i,j];
163 printf '\n';
164
165 printf 'Constraints\n';
166 printf 'Period Harvest Sign \n';
167 for {k in 1..K_PERIOD} {
168 printf '%5s %10.2f >= %.3f\n', k, sum {g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} Yield_Table_Vol[g,i,j,k] * XForestLand[g,i,j], Harvest_Min_Vol_Period;
169 }
170
171 # xbase (dbf) output
172 table Harvest{k in 1..K_PERIOD} OUT "xBASE" "HarvestArea1.dbf" "N(5)N(15,2)" : k ~ Period, (sum {g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} Yield_Table_Vol[g,i,j,k] * XForestLand[g,i,j]) ~ H_Area;
173
174 # xbase (dbf) read
175 set S, dimen 2;
176 table tab2 IN "xBASE" "HarvestArea1.dbf": S <- [Period, H_Area];
177 display S;
178
179
180
181
182 printf '\n';
183 printf 'Constraint\n';
184 printf 'Harvest Period\n';
185 printf 'Type AgeClass PrescMgt Period Value\n';
186 printf{g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR, k in 1..K_PERIOD}:'%5s %11s %11s %5s %10.2f\n', g,i,j, k, (Yield_Table_Vol[g,i,j,k] * XForestLand[g,i,j]);
187
188
189 printf 'Even_Flow_Constaint_Alpha (1-Alpha)\n';
190 printf 'Period Sign \n';
191 for {k in 6..K_PERIOD-1} {
192 printf "%s %10.2f <= %s\n", k, ((1 - Alpha) * sum {g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} Yield_Table_Vol[g,i,j,k] * XForestLand[g,i,j] - sum {g in G_STAND_TYPE,i in I_CULTURAL_PRES, j in J_MGT_YEAR} Yield_Table_Vol[g,i,j,k+1] * XForestLand[g,i,j]),0;
193 }
194 printf '\n';
195
196
197 # Forest / Land Constraints
198 printf '\n';
199 printf 'Total Area Constraint\n';
200 printf 'Type AgeClass PrescMgt Value Sign Total_Area \n';
201 printf '%5s <= %.3f\n',sum {g in G_STAND_TYPE, i in I_CULTURAL_PRES, j in J_MGT_YEAR} XForestLand[g,i,j], Total_Area;
202
203 printf 'Area\n';
204 printf 'Area Value Sign Areas_Stand\n';
205 for {g in G_STAND_TYPE} {
206 printf '%5s %10.2f <= %.3f\n', g, sum {i in I_CULTURAL_PRES,j in J_MGT_YEAR} XForestLand[g,i,j], Area_Stand_Type[g];
207 }
208
209
210 #DATA SECTION
211
212 data;
213
214 # Most of the data has been moved to dbf format
215
216 param MGT:=31;
217
218 param K_PERIOD:= 7;
219
220 param Alpha:= 0.20;
221 param Beta:= 0.20;
222
223 param Harvest_Min_Vol_Period:= 12000;
224
225 end;
226