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