LINEAR OPTIMIZATION OF FOREST MANAGEMENT FOR DYNAMIC RECURSIVE MODEL

O . T u r k o v s k a Postgraduate student* Research Assistant** E-mail: turkovska@gmail.com M . G u s t i PhD* Researcher** E-mail: kgusti@yahoo.com P . L a u r i PhD, Researcher** E-mail: lauri@iiasa.ac.at N . F o r s e l l PhD, Researcher** E-mail: forsell@iiasa.ac.at P . H a v l í k PhD, Senior Researcher** E-mail: havlik@iiasa.ac.at M . O b e r s t e i n e r PhD, Program Director** E-mail: oberstei@iiasa.ac.at *Department of International Information Lviv Polytechnic National University **Ecosystems Services and Management Program International Institute for Applied Systems Analysis Schlossplatz, 1, Laxenburg, Austria, A-2361 У дослідженні представлено розроблення моделі оптимізації лісокористування на основі лінійного програмування, яка може бути легко інтегрована у складні динамічні рекурсивні моделі і яка містить інструменти, що забезпечують врахування майбутнього стану лісу під час поточної заготівлі деревини за наявності рекурсивних обмежень. Представлено загальну структуру і алгоритм моделювання. Проведено порівняння результатів двох лісових моделей та визначено основні шляхи подальшого вдосконалення розробленої моделі лісокористування Ключові слова: лісокористування, лінійне програмування, рубки, вікова структура, вартість лісу, модель лісокористування, динамічна рекурсивна модель


Introduction
A forest management model is a combination of an economic model of decision making be forest owner and an ecological model of forest growth.Forest management models simulate the forest owner behavior, potential growth of forest and forest response to different treatments over time by using mathematical programming, heuristic methods and other simulation techniques.
Linear programming (LP) is the most used method and it has a long tradition in forest sector modeling [1].LP can solve efficiently complex large-scale optimization problems, which might be computationally intractable with other methods.
Forest management modeling can be characterized as an optimization problem because forest owner is looking for the best solution concerning forest planting, providing a favorable condition for forest growing and harvesting.The management technique depends on a goal which forest owner is willing to fulfill.Generally the goal is a trade-off between intensification of forest harvesting and carbon sequestration.
Management of forest differs significantly from management of agriculture.In general agriculture has one year managing cycle.The yield is harvested and then planted for the next year.Forest management requires more challenging planning.Harvesting regime has to take into account future yield of forest and environmental needs in order to keep forest under sustainable development over time.
Generally LP is applied for solving forest management when quantity length of management periods are known and optimization is intertemporal [1].This approach is not applicable for dynamic recursive models where optimization is under one period.Large scale economic models which include forestry as one of the economy sectors are commonly built based on the dynamic recursive scheme.However, such an approach has limitations when it comes to the optimization of the forest management as the decisions as of one period has a strong correlation with future  O. Turkovska, M. Gusti, P. Lauri, N. Forsell, P. Havlнk, M. Obersteiner, 2015 Математика и кибернетика -прикладные аспекты periods and decisions.In this paper we develop a new forest management modeling approach for recursive dynamic models that is taking into account current and future occurrences.

Analysis of published data and problem statement
Ecological models of forest growth can be divided in two classes.Traditional empirical (or statistical) growth models relay on growth functions and yield tables, which are estimated from historical forest inventory data [2].In practice, good forest inventory data is not always available, because collecting forest inventory data is time consuming and expensive.An alternative approach to model forest dynamics is process-based (or mechanistic) growth models, where woody biomass growth is modeled explicitly depending on the physical process such as photosynthesis, respiration, radiation and CO 2 concentration [3].The advantage of process-based growth models is that they are based on net primary productivity (NPP) estimates, which are currently easily available for all parts of the world.The disadvantage of process-based models is that they are often not as accurate as empirical models if we compare them to historical forest inventory data.
Spatial scale of forest management models ranges from stand or even tree level to landscape level.Stand level forest management models solve the optimal time to cut a single stand of trees.Usually it is assumed that the stand consist of even-aged trees (optimal rotation models), but there are also models that consider heterogeneous stands with mixed-species and uneven-aged trees (single tree models).Landscape-level forest management models determine how a given land area should be allocated between multiple stands.A classical solution to manage multiple stands is a normal forest [4].A normal forest is a steady state, which has a uniform distribution of homogenous stands with different age-classes.Each period the oldest stand is harvested and replanted such that it becomes the youngest stand.The normal forest is a desirable state of forest management, because it leads to constant sustainable yields over time.
Forest management models can be non-spatial or spatially explicit.Non-spatial models do not consider the specific location of forest management units, i.e., they do not utilize geographical information system (GIS) data.For example, they do not consider the transport distances between harvested forests and saw mills or energy plants.A spatially explicit forest management model includes the geographical dimension in the analysis.Spatially explicit forest management models have become more common in last 20 years with increased availability of GIS data and higher computing capacity of computers [5].
Spatially explicit landscape-level forest management models have typically a large number of area and ageclass data, i. e., the large-scale models.Because the size of dynamic optimization problems grows geometrically with the state space, this tends to cause a serious problem for the model solving.Usual way to overcome the state space problem is to use recursive optimization instead of intertemporal optimization.Recursive optimization sets some restrictions on the modeling, because in the recursive model decision makers are myopic and they do not anticipate the future states of the world correctly [6].However, recursive optimization allows to keep the solving times of the model in reasonable limits.

Purpose and objectives of the study
The aim of our study is development of a forest management model based on LP that can be easily integrated into a large-scale dynamic recursive model and contain the instruments which provide future consideration for harvesting plans under recursive limitations.
In accordance with the set goal the following research objectives are identified: -developing a forest management model that can provide temporal allocation of forest harvesting; -analysis of received modeling results and comparison of projected age structure with Global Forest Model (G4M); -identifying the possible future steps in order to improve adequacy of the developed forest management algorithm.

Description of forest management model FesT
We developed a spatially explicit forest management model (FesT) based on linear programming.The goal of the model is to simulate processes of forest management with the possibility of integration into a large-scale dynamic recursive model.FesT maximizes total benefit of forest value over simulation periods.The model operates on a regular geographic grid of 0.5x0.5 deg.It searches optimal time (period) of harvesting and optimal place (cell) of harvesting when producing amount of wood set externally.Currently the model runs 5 simulation periods i. e., till 2050.Inputs, parameters, variables and outputs of the model are shown in Table 1.Data and processes of the model are introduced on three scale levels i. e., country, grid cell and age class.Data which do not require a high spatial resolution or cannot be detailed on another scale are introduced on a country level, e. g. wood demand, planting and harvesting costs and wood prices.We can consider impact of natural conditions on forest growth using yield tables specified for every grid cell.FesT structure consists of three operational blocks: pre-calculation, optimization and updating (Fig. 1).The first block -pre-calculation, runs only once.It computes parameters on grid cell or age class levels for all simulation periods: rotation time, land expectation value, forest value, final cut area etc.Some parameters are computed only for first simulation period, in particular final cut area and will change for the next simulation periods.
In the second block the model optimizes harvested forest area in each age class of each cell where final cut area >0 in order to fulfill domestic wood demand of a country and maximize the total benefit of forest value.
In the third block the model updates the area of the age classes and the final cut area in order to prepare the input data for the next simulation period.Wood demand, forest value, benefit losses and delay costs are updated to their values at the next simulation period which were computed in the pre-calculation block.While forest area and final cut area in the age class are modified considering harvested area in the running simulation period.After this step the next simulation period starts.
Formally optimization problem of the model is expressed as following: ( ) ( ) where t -index which refers to simulation period (t=1…5 Harvested area in the grid cell and certain age class cannot exceed final cut area of those grid cell and age class.Harvested amount of wood on the country level must be equal or higher then demand but it cannot exceed demand more than by 1 %.
The forest value conception is useful for determining a harvesting time of the forest stand [7].Comparison of the forest value in the current time period and at some point in the future can help to take a decision when certain forest stand will be better to harvest.This approach is useful for the model because it can prevent harvesting at the current simulation period in order to perform it when this will be more profitable.According to the approach prices and costs remain the same over time.Considering that we can estimate forest value in the model for every simulation period before optimization starts.The forest value computation in the model is based on [7] and formally expressed as: (2 where T -time when forest in certain age class of certain grid cell will be harvested; wp -wood price in the country ($/m 3 ); gs c,i -growing stock of forest at certain age according to the yield table of certain grid cell (m 3 /ha); hc -harvesting costs in the country ($/ha); pc -planting costs in the country ($/ha); lev c -land expectation value at certain grid cell ($/ha); r -discount rate.
Land expectation value is used to estimate a value of bare land that may be used for forestry and find rotation time which will maximize the benefits periodic perpetual harvesting [7].In the model it indicates a forest management practice which will bring the highest benefits under current conditions when existing forest will be harvested: where rt c -rotation time in the grid cell c.We consider each age class as a separate forest stand.We assume that first simulation period is initial time point when owner of the forest stand is taking a decision when to harvest the forest stand.The forest stand can be harvested in current simulation period or harvest can be postponed to the next simulation periods.Consequently T=0 for the first simulation period and increases by 10 years for each follow-Fig.1.General structure of FesT ing simulation period.Knowing the forest value of current simulation period is not enough to evaluate more profitable harvesting time.We determine two parameters -benefit losses and delay costs which are the drivers of harvesting decision.
Benefit losses (bl t,c,i ) determines a benefit from forest harvesting which forest owner loses if harvesting occurs in the running simulation period (t): If forest value in running simulation period (t) is higher than in the following simulation period (t+1) it means the forest owner will not lose any profit if he harvests now.If the forest value is higher in the following period then we assume that the forest owner will lose the profit which he could earn in the future.
Delay costs (dc t,c,i ) deal with other side of the problem.It determines how much the forest owner will lose if he does not harvest in the running period (t): If forest value in running simulation period (t) is lower than in the following simulation period (t+1) it means the forest owner will not lose any profit if he postpones harvest.If the forest value is higher in the following period then we assume that the forest owner will lose the profit if he postpones harvesting because in the following period income will be lower than in the running period.

Results of developed forest management model
Described approach of forest management modeling can be used globally.We test and analyze model's behavior on one country -Ukraine.Current version of the model is focusing on forest management activities and does not include detailed forest growth modeling.FesT is linked with Global Forest Model (G4M) which provides the yield tables for Ukrainian forest, area of the forest and its age structure distribution [8].FesT is also linked with Global Biosphere Management Model (GLOBIOM) which provides data on wood price, harvesting and planting costs [9].Historical data for wood demand are based on statistics of Food and Agriculture Organization (FAO) [10].Wood demand projections until 2050 keep a dynamic of the historical data.
We run the model for five simulation periods to test how benefit losses and delay costs drive spatial and temporal distribution of harvesting.At the (Fig. 2-6) is shown comparison of forest area distribution among age classes on the country level during all simulation periods: The forest area is represented by three groups i. e., total forest area, final cut area and harvested area.Total forest area contains all forest in the country distinguished among age classes before optimization would occur in the simulation period.The final cut area contains forest which can be harvested in the country in the simulation period.The harvested area contains forest which was harvested from the age class in the simulation period.The final cut area is estimated due to rotation time.The forest younger than rotation time cannot be harvested.In the model we calculate rotation time to maximize annual increment which prevent forest from too intensive felling.At the (Fig. 2-6) we can observe how forest area is changing through time between age classes.We assume that at the beginning of the simulation period we plant forest an area equal to the total harvested area in the previous period.If harvesting does not occur in the run-Восточно-Европейский журнал передовых технологий ISSN 1729-3774 5/4 ( 77 ) 2015 ning simulation period at a certain age class then forest shifts to the next age class in the following simulation period with the same area.This case can be observed on the young forest transition.If harvesting occurs then forest shifts to the next age class taking into account the harvested area.Final cut area for the following simulation period is estimated after the forest area changes due to harvesting in the running period.

Analysis of the FesT results
Current version of the model constraints harvest exceeding domestic demand even if harvesting is profitable.Computation of forest value (2) shows that it increases until forest becomes a middle age and after this point forest value starts decreasing.Such forest value dynamics occurs because we consider growing stock but do not consider share of saw logs.It means that benefit losses are high for the young forest and delay costs are high for mature forest.As a result we see that it is more profitable to harvest forest now than wait.Every simulation period the model is harvests the forest in order to reduce the delay cost, consequently it harvests mature and overmature forest.This harvesting practice differs from is the harvesting regime used in Ukraine.Most final cuts occur in 40-80 years and the number of overmature forest is increasing [12].It is happening because generally the quality of wood is better in this age range.To consider this case the method of calculation of forest value must be improved by taking into account additional characteristics of trees, e. g. diameter.Including data about group species (coniferous, broadleaves) and estimating rotation time for each group will also remove unnecessary high constraints on harvesting age.
Rotation time is an important tool to control harvesting volumes and keep stable growing stock increment for the future periods.In reality rotation time varies for different species and natural conditions.Computation of rotation time at least for different groups of species will provide more accurate limitations on harvesting intensity.
To validate performance of FesT, in particular the dynamics of age-class distribution we compare the model results with G4M We run both models from 2000 till 2050 with the same inputs and compare the projected age-class distribution averaged over the country after the last simulation period (Fig. 7).
G4M is a spatially explicit global forest model which simulates forest management and land-use changes over time and their response to climate policies [8].Both models are spatially explicit, simulate dynamics of age-class distribution and harvesting levels are driven by domestic wood demand.However, G4M is not an optimization model.Methods of forest management simulation differ significantly in G4M and FesT.G4M sorts grid cells which contains forest by mean annual increment (MAI), biomass, forest area and population density.The model assumes that high-productive large forests that are located in the grid cells with high population density (closer to markets) are more profitable for harvesting.G4M selects management type for a forest in each grid cell taking into account its mean annual increment and net present value.For each grid cell G4M computes three types of rotation time: maximizing mean annual increment, keeping current biomass and maximizing biomass.Next step is selection of rotation time for the forest according to wood demand.The rotation time is selected in the range from minimum rotation time providing sustainable harvest (maximizing mean annual increment) and maximum rotation time when forest reaches maximum biomass.G4M provides forest management in such a way to bring the forest to a "normal" state (equal area of all age classes) after one rotation period.The model harvests in the forest with higher MAI and compares harvesting level with wood demand.If the harvest does not match the demand the model adjusts rotation time and uses forest with lower MAI (if it is necessary) to get the difference between demand and supply lower than ±1 % [8,11].This approach prevents G4M from using all forest area for wood production, if wood demand can be satisfied by harvesting only forest with higher MAI, and remains forest with lower MAI unused.The algorithm of forest stands transition from one ageclass to another is similar in both models that explain why the models finish the simulation with very close age-class distributions.However, there is some divergence that has to be analyzed.
At first six age classes area is larger in FesT's projection.G4M and FesT apply the same rule for forest planting: the area of planted forest in the following simulation period equals to the harvested area in the running simulation period.Wood demand is the same for both models so they have to harvest the same amount of wood.But the models are free to select where exactly this wood comes from.
The age classes starting from the 7th are the ones where harvest occurs in both models.In these age classes G4M has larger forest area than FesT.Difference in harvest modeling explains the difference in forest area.In G4M harvesting is performed in cells with forests with higher MAI leaving some cells unused that results in larger forest area in older age classes when averaged over the country.The significant difference in forest area occurs in the last age class.As it was mentioned before G4M harvests the forests with higher MAI at first then if necessary it shifts to forests with lower MAI.Therefore, forests with low MAI continue to grow and could stay untouched during all simulation periods.This causes a large forest area in the oldest age class.In FesT forest in the oldest age class is always harvested due to very high delay costs.
This comparison clearly shows that applying financial mechanisms to regulate forest management considers the decision on time and location of harvest from one side.The current version of FesT does not contain explicitly drivers which can direct the model to account for forest productivity and distance to markets.Increasing of harvested area is not always the best solution for fulfilling the demand.Therefore, these drivers must be included into the FesT to make its decisions more reliable and balanced.Natural mortality is necessary to integrate for better representation for natural processes of forest and more accurate estimation of forest available harvesting.
Harvested area in FesT is constrained by wood demand, harvesting and planting costs.The impact of harvesting and planting costs is not evident.The reason for that is a high growing stock which compensates the costs.Therefore, the model does not minimize harvested forest area.
In FesT share of harvested forest area in an age class depends on forest value, benefit losses and delay costs (both are negative values) (Fig. 8).The delay costs are higher for older forest because future increment does not compensate discounted profit.FesT harvests in those age classes where the delay cost are the highest.
Forest value gets its highest point at the same age as growing stock.Delay costs are driven by current increment of growing stock.With age current increment decreases that cause increasing of respective delay costs.In the age of 170 growing stock reaches its highest value and starts to decrease.At the same time current increment after 170 years becomes negative.Consequently, delay costs are the highest at 170 years and after they decrease.Through comparison of the forest values we can deal with recursive model limitations and evaluate better time for final cutting the forest.This approach expresses financial consequences of final cut decision in the running simulation period.It does not estimate the harvesting time precisely but it directs the model to perform it in better period, i. e. age class.However, a few improvements can be made in the forest value calculations.One of them is including qualitative characteristics of the harvestable wood.The results of the model show that harvesting occurs in old forest.The reason for it is avoiding high delay costs, but it does not take into account quality of the old wood.Consequently, consideration of wood quality in different ages or different price will make the harvesting strategy of the model more realistic.

Conclusions
We introduced a new approach for modeling forest management under recursive dynamic model.Concepts of forest value, harvest losses and delay costs provide harvesting decision which considers future state of the forests.The forest growth modeling is introduced in the FesT simplistically.We applied a generalized yield tables for all Ukrainian forests without distinguishing among tree species.Influence of natural conditions on forest productivity is taken into account by applying geographically explicit approach.In the model forest growth and forest management modeling are split between different processing blocks.Therefore, the optimization process is not overloaded by processing of additional data and forest growth modeling can be improved without interference into optimization.
FesT yields the same dynamics of forest area in the age classes over time as G4M.Different approaches of forest management modeling in the models explain the divergence between area distributions over age classes.
However, there are several aspects in the model that can be improved.Currently benefit losses do not impact harvesting decision.Distinguishing two categories of wood: saw logs and rest wood will force the model to wait longer before harvesting in order to get enough saw logs to fulfil demand.Therefore, this will increase the benefit losses and their influence on harvesting decision.Distinguishing tree species according to growth rate will provide a different dynamics for fast and slow growing tree species and will allow estimation of the values of rotation time, benefit losses and delay costs more accurate.Natural mortality of forest is introduced in a simple way and only in the oldest age class.This process must be introduced in all age classes but with different intensity.Implementation of the abovementioned improvements will make modeling of forest management closer to reality.

Fig. 7 .
Fig. 7. Projection of age structure distribution of Ukrainian forests in 2050 by two models

Table 1
Inputs, parameters, variables and outputs in the FesT Rotation time (time passing from planting of a forest stand till it is harvested) is computed considering productivity of the forest of certain grid cell.Forest area in each grid cell is distributed among forest stands of different age grouped, e. g. 10 years (age classes).The age classes where harvesting can occur are defined through rotation time, area of the age classes determines allowed final cut area.In the current version of the model age class width is fixed to 10 years and there 22 age classes.The first age class contains forest from just planted up to 10 years old and the last one contains forest 210 years and older.
age class under a certain grid cell (i=1… 22); H t,c,i -harvested area (ha) in age class i of grid cell c in the running simulation period t; fv t,c,i -value ($/ha) of 1 ha of forest stand in age class i of grid cell c in the running simulation period; bl t,c,i -benefit losses ($/ha) in the age class i of grid cell c if the forest stand is harvested in the running simulation period; fc t,c,i -area in age class i of grid cell c where final cuts are allowed in running simulation period; dc t,c,i -delay costs in the age class i of grid cell c if the forest stand is harvested in the next simulation period.
); c -index which refers to a certain grid cell; i -index which refers to an