Quality of Threshold Accepting Heuristic Search Results WhenAttempting to Solve Forest Harvest Scheduling Problems

Abstract

Threshold accepting is a s-metaheuristic, local search process that moves from one solution to another within the feasible region of the solution space via a random change to one or more elements of the solution. Threshold accepting can be further characterized as an aspirational combinatorial optimization process that does not guarantee optimality. The quality of outcomes from a threshold accepting search process varies when applied to forest harvest scheduling problems depending on the parameter value assumptions and sub-processes employed. Two relatively small but realistic case study forests are subjected to four management scenarios and outcomes are examined to illustrate how the quality of solutions may differ when the parameter values and processes employed within threshold accepting are adjusted. Statistically significant improvements in solution quality were generally evident with a slowing of the rate of change in the threshold value and the enhancement of the search process by using 2-opt moves and search reversion. While it was rarely observed, the threshold accepting heuristic search process located the optimal solution of most of the problems modeled. In cases where the problems involved maximizing an economic objective, about 47% of the heuristic search solutions had an objective function value that was within 1% of the optimal solutions.

INTRODUCTION

Many contemporary forest harvest scheduling problems require the use of discrete decision variables and mixed integer problem formulations. The desire to control the timing, placement, and size of final harvest activities is one example. The development of feasible and efficient solutions for mixed integer programming problems can be arduous with respect to the computing time required. Mathematical programming using exact or algorithmic methods represents the act of locating an optimal solution to a problem which was expressed using mathematical notation (Killen, 2021) using procedures that can guarantee locating the optimal solution (Romanycia and Pelletier, 1985). These methods are the preferred approach for forest harvest scheduling efforts, however, the challenges related to environmental concerns have prompted the need for integer decision variables, and thus some contemporary forest harvest scheduling problems are difficult or impossible to solve with exact methods (Brodie and Sessions, 1991). Due to these challenges, heuristic search processes have been suggested as potential methods to assist with forest harvest scheduling efforts.

S-metaheuristic search processes operate one level above basic move selection search processes (e.g., hill-climbing or random search) in navigating through a solution space (Albus, 1981). These types of search processes capitalize on insight or knowledge gained during recent search activity (Romanycia and Pelletier, 1985). Originally proposed by Dueck and Scheuer (1990), a threshold accepting search process is considered a s-metaheuristic (Talbi, 2009), a unidirectional approximation algorithm for solving complex problems (Siedentopf, 1995), and a refined local search algorithm (Winkler and Fang, 1997). One solution to a problem is maintained in computer memory at any one point in time during the search process, and an adjustment (a move within the neighborhood of the solution space) is made to convert this solution to a slightly different solution in the local neighborhood of the solution space. As with simulated annealing, a potential move is randomly selected from the set of possible adjustments to the current solution. In forest harvest scheduling, for example, a move might consist of randomly selecting a stand of trees, then randomly selecting a different year in which these trees will be harvested. The objective function value and the constraints are then assessed, and the move is potentially accepted.

Threshold accepting is very similar in concept and structure to simulated annealing, except for the manner in which it decides whether inferior moves suggested by the search process are acceptable. As others have noted (Nissen and Paul, 1995), a threshold accepting search process is a simplification of a simulated annealing search process. Both search processes attempt to change a current feasible solution to a new feasible solution through a random perturbation of the current solution. When the proposed change leads to a higher quality feasible solution (as reflected in an improvement in the objective function value), both search processes automatically accept the proposed change, and the new solution becomes the current solution (Fig. 1). When the proposed change leads to a lower quality feasible solution, the search processes differ in how they decide whether to accept the change (even though it may lead to a lower quality solution) or to reject the change. When using a simulated annealing search process a change to a lower quality solution is acceptable with a given probability that is a function of how long the search has been conducted (as reflected by the current temperature of the annealing process) and the difference in objective function values between the proposed feasible solution and the current (or often the best) feasible solution. The temperature value within a simulated annealing search process begins relatively high and decreases as the search progresses. The probability of acceptance then decreases as the temperature value decreases (as the search progresses) and as the difference in objective function values increases.

Fig. 1.

A general flow chart of the threshold accepting heuristic search process.

The simplification employed within threshold accepting involves replacing the probability of acceptance with a straightforward assessment of the difference between the objective function values of the proposed, temporary solution and the current (or best) solution. So, if a proposed feasible solution has an objective function value that is lower in quality than the objective function value of the current (or best) solution that is held in memory, yet the difference between the objective function values of these two solutions (the proposed solution and either the current or best solution) is smaller than the current threshold value, the temporary solution is deemed acceptable, and it becomes the current solution. However, if the difference between the objective function values of the two solutions (the proposed solution and either the current or best solution) is larger than the current threshold value, the temporary random adjustment is discarded (the proposed solution is rejected) and the search reverts back to the previous feasible solution to the problem. The threshold value within a threshold accepting search process begins relatively high and decreases as the search progresses. At some point the threshold becomes zero (0) and the search process becomes a local hill-climbing search. Threshold accepting may be conceptually easier to comprehend than simulated annealing because the threshold value is represented with the same units as the objective function value of the problem being solved. For example, if an objective function is designed to assess the value of a harvest schedule in Japanese Yen, the threshold value is also represented as Japanese Yen rather than a probability of acceptance.

From the description provided here it is evident that threshold accepting is a somewhat discriminatory, somewhat forgiving, stochastic hill-climbing search process. In part due to its speed and simplicity, threshold accepting has been used in forestry to address complex forest harvest scheduling problems (Bettinger et al., 2002; Bettinger et al., 2003; Bettinger and Boston, 2008; Zhu and Bettinger, 2008; Li et al., 2010; Bettinger et al., 2015; Akbulut et al., 2017) and has been used to address complex problems in a wide variety of other domains of commerce (e.g., Tarantilis et al., 2004; Perea et al., 2008). Other tangential issues in forest management can be addressed as well with a heuristic approach such as this. For example, with the assistance of a threshold accepting search process, Marshall et al. (2003) illustrated how smaller land units constructed within a geographic information system allow greater flexibility in the scheduling of management activities with a threshold accepting search process. And, Coulter et al. (2006) illustrated how a threshold accepting search process could be used to maximize the value of a schedule of forest road maintenance activities that are subject to budgetary constraints.

The objective of this paper is to illustrate how the outcomes from a threshold accepting search process may change as the parameters of the search change, and as sub-processes are added to facilitate diversification and intensification of the search. Statistical tests are employed to assess whether two sets of solutions developed utilizing different search assumptions produce solutions that are significantly different. The general hypothesis is that the quality of solutions generated from a threshold accepting search process is the same regardless of differences in how the search is conducted.

METHODS

Two hypothetical forests are used to illustrate the differences in outcomes when threshold accepting parameters are adjusted. Both forests are subject to problems where there is a maximization (net discounted revenue) and a minimization (deviations from a target wood flow value) objective. The Lincoln Tract (Bettinger et al., 2017) is a hypothetical coniferous forest located in the western United States. It is composed of 87 stands covering 1,841.5 hectares, although three stands are unavailable for harvest as they represent predominantly riparian areas. The remaining 84 stands encompass 1,788.4 ha. Douglas-fir (Pseudotsuga menziesii) is the most dominant tree species within the forest. The time horizon is 30 years, and each of the 6 time periods are 5 years long. Harvests are assumed to occur, for modeling purposes, in the middle of each time period. The volumes are represented in thousand board feet (MBF), a traditional unit of sawtimber in the United States. A harvest target of 13,950 MBF per time period is assumed in the minimization problem, an amount that was estimated based on the Hanzlik formula (Hanzlik, 1922). One management action (final harvest) is modeled, and the minimum final harvest age was assumed to be 35 years. The maximum size of a final harvest area is assumed to be 48.6 ha (legal maximum in Oregon). The green-up period after a final harvest has been scheduled is assumed to be 5 years (one time period).

The Jones Tract is a hypothetical coniferous forest located in the southern United States. It is composed of 81 stands covering 1,053.1 hectares, although 16 stands are unavailable for harvest as they represent predominantly wetland areas. Thus 65 stands of trees are available for harvest activities (867.4 ha). Loblolly pine (Pinus taeda) is the most dominant tree species within the forest. The time horizon is 20 years, and the time periods are each 5 years long. Harvests are assumed to occur, for modeling purposes, in the middle of each time period. The volumes are represented in tons (907.2 kg), a traditional unit of weight in the United States. A harvest target of 19,000 tons per time period is assumed in the minimization problem, an amount estimated based on the Meyer amortization method (Meyer, 1952). One management action (final harvest) is modeled, and the minimum final harvest age was assumed to be 22 years. The maximum size of final harvest area is assumed to be 48.6 ha. The green-up period after a final harvest has been scheduled is assumed to be one time period.

Four forest management problems are addressed in this work:

• 1. Maximize discounted net revenue from final harvests and employ unit restriction harvest adjacency (URM) constraints (Murray, 1999) to control the timing and placement of final harvests

• 2. Maximize discounted net revenue from final harvests and employ area restriction harvest adjacency (ARM) constraints to control the timing and placement of final harvests

• 3. Minimize squared deviations between the scheduled harvest volumes in each time period of the time horizon and a target harvest volume and employ a unit restriction harvest adjacency model to control the timing and placement of final harvests

• 4. Minimize squared deviations between the scheduled harvest volumes in each time period of the time horizon and a target harvest volume and employ an area restriction harvest adjacency model to control the timing and placement of final harvests

Minimizing the squared deviations between scheduled harvest volumes and a target harvest volume leads to a more even scheduled harvest volume across all time periods than simply minimizing the straight differences between scheduled harvest volumes and a target harvest volume. The two objective function values for the management problem formulations are:

  

$Maximize\sum _{t=1}^T\sum _{n=1}^N\left(\frac{x_{tn}{r}_{tn}}{p}\right)$(1)

  

$Minimize\sum _{t=1}^T{\left(H_t-TV\right)}^2$(2)

These problems include an accounting row or function

  

$\sum _{n=1}^N\left(x_{tn}{v}_{tn}{a}_n\right)-H_t=0\forall t$(3)

And only one (at most) final harvest can be assigned to each stand during the time horizon.

  

$\sum _{t=1}^T{x}_{tn}1n,x_{tn}\left\{0,1\right\}$(4)

When the unit restriction model of adjacency is assumed, the following constraint applies when a final harvest is applied to a stand (n) during time period t:

  

$x_{tn}+\sum _{z=lw}^{uw}\sum _{m{O}_n}^M{x}_{zm}1n,t$(5)

  

$lw=t-\left(g-1\right)$(6)

  

$iflw<1,lw=1$(7)

  

$uw=t+\left(g-1\right)$(8)

  

$ifuw>T,uw=T$(9)

Here, variable z represents the length of time known as the green-up period. In this case, each stand (n) and each neighboring stand (m) are assessed as a pair, and at most only one of these is scheduled a potential final harvest within the green-up window defined by time period t. When the area restriction model of adjacency is assumed, the following constraint applies:

  

$x_{tn}{a}_n+\sum _{z=lw}^{uw}\sum _{m{O}_n{S}_n}^M{x}_{zm}{a}_{m}MA$(10)

In this case, each stand (m) that is adjacent to the focal stand (n), along with any other management units adjacent to m, and their neighbors, and so on, are assessed as a sprawling cluster of potential final harvest areas all scheduled within the green-up window defined by time period t.

For the maximization problems, the following wood flow constraints apply:

  

$\left(\left(\sum _{t=1}^T{H}_t\right)/T\right)-AVGV=0$(11)

  

$H_t\left(1-\right)AVGVt$(12)

  

$H_t\left(1+\right)AVGVt$(13)

where

a_n = area of stand n

AVGV = average scheduled volume across all time periods in the time horizon

β = allowable deviation (decimal percentage, where 25% = 0.25) from average scheduled wood volume

g = green-up window, expressed as number of time periods t

H_t = scheduled harvest volume during time period t

i = interest rate for discounting purposes

lw = the lower limit (time periods) on the period of time within a green-up window of a proposed final harvest

m = a neighboring stand (management unit) to stand n

M = the total number of neighboring stands to stand n

MA = maximum final harvest area size

n = a single timber stand from the model forest

N = the total number of timber stands in the model forest

O_n = the set of neighboring stands to stand n

  

$p={\left(1+i\right)}^{\left(t-1\right)\text{*}5+\left(0.5\text{*}t\right)}$

r_tn = potential revenue for scheduling stand n for a final harvest during time period t

S_n = the set of all stands adjacent to stands in set O_n and all stands adjacent to neighbors of neighbors, and so on (Murray, 1999)

t = a single time period

T = the total number of time periods

TV = target harvest volume per time period

uw = the upper limit (time periods) on the period of time within a green-up window of a proposed final harvest

v_tn = volume per unit area available for harvest in stand n during time period t

x_tn = a binary decision variable indicating whether (1) or not (0) stand n is scheduled for a final harvest during time period t

z = a period of time within a green-up window of a proposed final harvest

The interest rate for discounting purposes (i) is assumed to be 5% in the maximization problem, and timber volumes per time period are allowed to vary (β) by 25% from the average scheduled volume per time period when scheduled.

The maximization problems were subjected to the branch and bound search process within Lingo version 21 (LINDO Systems Inc., 2024). The minimization problems were subjected to quadratic branch and bound search process within Lingo version 21. The threshold accepting heuristic search process, developed using Visual Basic 2012, was employed to generate solutions to these problems. Maximization problems required less than one second to solve with Lingo. Minimization problems, due to the desire to achieve exact even flow, were allowed to run for 50 to 110 hours before the quadratic branch and bound search process was interrupted. Even though the termination point for these searches was unpredictable, we assume the results from these are effective representations of the optimal solution. Heuristic search required about 1 second per instance (run) for the least diligent search and about 10 minutes for the most diligent search.

The general parameters of a threshold accepting search process include:

• • Initial threshold value

• • Rate of change of the threshold value

• • Number of successful moves (iterations of the model) to conduct before changing the threshold value

• • Number of unsuccessful moves to attempt before changing the threshold value

The latter of these parameters can be based on consecutive or non-consecutive unsuccessful moves since the last time that a successful change was made to the current solution during a search. In the work presented here, it is assumed that the number of consecutive unsuccessful attempts triggers a change in the threshold value. In any event, this element of threshold accepting (changing the threshold due to a large number of unsuccessful attempts to change the current solution) is necessary for the search process to eventually terminate. Unfortunately, this also requires extra coding to track the recognition of unsuccessful attempts (gray boxes in Fig. 1). Methods for utilizing the unsuccessful move strategy were described in Siedentopf (1995). Coulter et al. (2006) also used the unsuccessful count as a way of terminating the search process. The rate of change of the threshold value could either be static (i.e., change the threshold by X amount each time), adaptive (e.g., change the threshold value by a variable amount depending on the status of the search through the solution space), or non-linear (e.g., multiply the current threshold value by a rate of change (such as 0.99) to create a new threshold value). In this work, we are assuming that a non-linear rate of change will be employed each time the threshold value is changed. As one might gather, there are an infinite number of combinations of search parameters. Although Pukkala and Heinonen (2006) suggest that the parameters of a threshold accepting heuristic search process might be optimized for specific problems using a Hook and Jeeves direct search approach, here the parameters were selected specifically within certain ranges to illustrate the effect they have on the outcome of the search. And finally, some search processes allow infeasible solutions to a problem to be maintained by penalizing each constraint violation, but for the purposes of the work presented here, it is assumed this is not the course of action chosen.

With each of these four problems noted above, 36 different combinations of threshold accepting parameters are employed. These were selected to represent significant changes in the magnitude of parameter values. For example, for the Lincoln Tract, the initial threshold ranges from 1,000,000 to 5,000,000 US dollars, and the number of iterations per threshold is assumed to be either 1 or 100 (Table 1). For the minimization problems, the objective function value represents squared deviations of scheduled harvest volumes from target harvest volumes, and the initial threshold values vary from 1,000,000 to 5,000,000 units (squared deviations). Since the Jones Tract is smaller than the Lincoln Tract, when maximizing discounted net revenue, the largest value obtainable is around $2,000,000 (US dollars), therefore the initial threshold values when maximization processes are applied to the Jones Tract are 200,000, 350,000, and 500,000 (Table 2). As others have noted (Perea et al., 2008), threshold accepting requires calibration of a number of parameters whose values vary by problem instance. In the work presented here, it is illustrated why calibration might be necessary. Preliminary trials indicated that the initial threshold value for the maximization problems should be between 10 and 25% of the ultimate final objective function value, otherwise relatively small initial threshold values result in a search process much like hill-climbing, and relatively large initial threshold values result in a search process much like random search. In either case, there is a chance that the heuristic search will become stranded on poor quality local optima and be unable to escape to higher quality areas of the solution space.

Table 1.Taxonomy of parameter alternatives used in assessing the threshold accepting heuristic when applied to the Lincoln Tract (maximization and minimization problems) and the Jones Tract (minimization problem).

	Alternative
	1	2	3
Rate of change	0.99	0.999	0.9999
	(1)	(2)	(3)
Iterations per threshold	1	100
	(a)	(b)
Unsuccessful iterations per threshold	500	1000
	(α)	(β)
Initial threshold value	1,000,000	3,000,000	5,000,000
	()	(●)	()

Example: Alternative 1aα uses a rate of change of 0.99, 1 iteration per threshold, 500 unsuccessful iterations per threshold, and an initial threshold value of 1,000,000.

Table 2.Taxonomy of parameter alternatives used in assessing the threshold accepting heuristic when applied to the Jones Tract (maximization problem).

	Alternative
	1	2	3
Rate of change	0.99	0.999	0.9999
	(1)	(2)	(3)
Iterations per threshold	1	100
	(a)	(b)
Unsuccessful iterations per threshold	500	1000
	(α)	(β)
Initial threshold value	200,000	350,000	500,000
	()	(●)	()

Example: Alternative 1aα uses a rate of change of 0.99, 1 iteration per threshold, 500 unsuccessful iterations per threshold, and an initial threshold value of 200,000.

Each of the four problems noted above were applied to the two model forests, and for each of the 36 sets of threshold accepting parameters assumed, 200 independent runs of the threshold accepting heuristic search process were generated, each producing a feasible forest harvest schedule. In general, 1-opt moves are employed, where one stand is randomly and temporarily selected from the current feasible solution and a change to the harvest timing for that stand is also randomly and temporarily selected. The objective function value and feasibility are then assessed to determine whether the temporary change to the harvest timing of the stand is acceptable. Three enhancements are applied in conjunction with the set of parameter values representative of the most diligent search. The first involves the use of 2-opt moves (Bettinger et al., 1999; Caro et al., 2003), where the harvest timing of two randomly selected stands are temporarily exchanged, after which the objective function value and feasibility are then assessed to determine whether the temporary changes to the harvest timing of the two stands are acceptable. Methods for utilizing 2-opt moves within threshold accepting have been described previously (Nissen and Paul, 1995; Tarantilis and Kiranoudis, 2002; Tarantilis et al., 2004). Since 2-opt exchange moves interject no diversity into the solution (i.e., the number of final harvests per time period remains the same), 1-opt moves need to be intermixed into the search process. Here, we employ a repeating pattern of 50 1-opt moves followed by 10 2-opt moves. The second additional process involves the use of search reversion, which has been shown to be of value in forest harvest scheduling efforts (Bettinger et al., 2015). In search reversion, after a certain number of moves the current solution is replaced by the best solution stored in memory. Here, we assume that search reversion occurs every 50 moves. The third employs both 2-opt moves and search reversion. An additional 200 independent runs were generated each of these three enhancements for each forest / management problem. Thus, for each forest (2) and each management problem (4), 7,800 heuristic search solutions were generated. In total, 62,400 solutions (harvest schedules) were generated by the threshold accepting heuristic search process. These represent a small sample of the potential harvest schedules for the Lincoln Tract (7⁸⁴) and the Jones (5⁶⁵), where the population of harvest schedules, whether feasible or not, when using discrete decision variables is

  

$\text{Number of different solutions }={\text{ Options for each stand}}^{\text{Number of stands}}$

The options for each stand include the choice to employ a final harvest in any period of the time horizon and the choice to not employ a final harvest at all. Thus, there are 7 options for each stand in the Lincoln Tract even though there are only 6 time periods in the time horizon.

As entry into the random number list of the personal computer employed (12^th generation Intel^® Core™ i9-2900, 2.4 GHz processor and 64 GB RAM) was based on the computer's clock, each of the 200 runs within each set can be considered independent samples from the larger population of potential solutions to the forest management problem (Golden and Alt, 1979; Los and Lardinois, 1982). Statistical tests were therefore conducted in order to assess whether significant differences can be asserted when comparing two sets of objective function values. A set of objective function values is defined by a combination of the forest considered, the management problem applied, and the parameters / processes employed during the search process. The two-tailed Student's t-test was applied to assess whether statistically significant differences exist between pairs of sets of solutions, as the distribution of objective function values within each set was approximately normally distributed.

RESULTS

General Results

From the observations of objective function values generated by the heuristic search process, it is obvious that slowing the rate of change in the threshold improves solution values for both the maximization (Fig. 2-5) and the minimization (Fig. 6-9) problems. With each set of three results (e.g., 1aα, 2aα, 3aα) improvements in the median and interquartile range of solution values are evident (i.e., 3aα results are better than 2aα results, which are better than 1aα results). Further, in the far-right side of each figure, the three sets of solution values that represent enhancements upon the most diligent search (3bβ) generally suggest that these enhancements are necessary to produce a solution near the optimal solution to each problem. In the maximization problems for both the Lincoln and Jones Tracts, it is also evident that sets of solutions generated with 100 iterations per threshold value were generally of higher quality than those which utilized only one iteration per threshold (more clearly evident in the group beginning with 1bα and extending to 3bβ in Fig. 5). In addition, the interquartile range of these solutions is much smaller than others where only one iteration per threshold was employed. These findings are not as evident with the minimization problem results. Furthermore, when comparing URM and ARM outcomes (e.g., Fig. 4 and Fig. 5), the ARM outcomes have higher quality objective function values due to the added flexibility the ARM constraint allows in assigning a harvest period to a timber stand. Some other general outcomes include the following:

Fig. 2.

Box (interquartile range) and whisker (full range) graph illustrating the threshold accepting heuristic search results for the Lincoln Tract and forest management problem 1 (maximization of discounted net revenue and unit restriction adjacency of final harvests). 3bβx represents the most diligent search parameters along with 2-opt moves, 3bβy represents the most diligent search parameters along with search reversion, and 3bβz represents the most diligent search parameters along both 2-opt moves and search reversion.

Fig. 3.

Box and whisker graph illustrating the threshold accepting heuristic search results for the Lincoln Tract and forest management problem 2 (maximization of discounted net revenue and area restriction adjacency of final harvests). 3bβx represents the most diligent search parameters along with 2-opt moves, 3bβy represents the most diligent search parameters along with search reversion, and 3bβz represents the most diligent search parameters along both 2-opt moves and search reversion.

Fig. 4.

Box and whisker graph illustrating the threshold accepting heuristic search results for the Jones Tract and forest management problem 1 (maximization of discounted net revenue and unit restriction adjacency of final harvests). 3bβx represents the most diligent search parameters along with 2-opt moves, 3bβy represents the most diligent search parameters along with search reversion, and 3bβz represents the most diligent search parameters along both 2-opt moves and search reversion.

Fig. 5.

Box and whisker graph illustrating the threshold accepting heuristic search results for the Jones Tract and forest management problem 2 (maximization of discounted net revenue and area restriction adjacency of final harvests). 3bβx represents the most diligent search parameters along with 2-opt moves, 3bβy represents the most diligent search parameters along with search reversion, and 3bβz represents the most diligent search parameters along both 2-opt moves and search reversion.

Fig. 6.

Box and whisker graph illustrating the threshold accepting heuristic search results for the Lincoln Tract and forest management problem 3 (minimization of discounted net revenue and unit restriction adjacency of final harvests). 3bβx represents the most diligent search parameters along with 2-opt moves, 3bβy represents the most diligent search parameters along with search reversion, and 3bβz represents the most diligent search parameters along both 2-opt moves and search reversion.

Fig. 7.

Box and whisker graph illustrating the threshold accepting heuristic search results for the Lincoln Tract and forest management problem 4 (minimization of discounted net revenue and area restriction adjacency of final harvests). 3bβx represents the most diligent search parameters along with 2-opt moves, 3bβy represents the most diligent search parameters along with search reversion, and 3bβz represents the most diligent search parameters along both 2-opt moves and search reversion.

Fig. 8.

Box and whisker graph illustrating the threshold accepting heuristic search results for the Jones Tract and forest management problem 3 (minimization of discounted net revenue and unit restriction adjacency of final harvests). 3bβx represents the most diligent search parameters along with 2-opt moves, 3bβy represents the most diligent search parameters along with search reversion, and 3bβz represents the most diligent search parameters along both 2-opt moves and search reversion.

Fig. 9.

Box and whisker graph illustrating the threshold accepting heuristic search results for the Jones Tract and forest management problem 4 (minimization of discounted net revenue and area restriction adjacency of final harvests). 3bβx represents the most diligent search parameters along with 2-opt moves, 3bβy represents the most diligent search parameters along with search reversion, and 3bβz represents the most diligent search parameters along both 2-opt moves and search reversion.

• • With respect to maximization problems (1 and 2) when applied to both case study forests, a relatively large range of final solution values were generally observed when the rate of change in the threshold value was relatively fast (0.99) as compared to slower rates of change.

• • With respect to both the minimization and maximization problems (1-4) when applied to both case study forests, repeatable high quality solutions, as reflected in a narrow range of solution values and a relatively short interquartile range, a slow rate of change in the threshold value seems necessary.

• • With respect to the minimization problems (3 and 4) when applied to both case study forests, threshold accepting was generally not able to perform well with the search parameters employed without the addition of enhancements (2-opt, reversion, or both).

Increased diligence of the threshold accepting search process is achieved through (a) a higher initial threshold value, (b) a slower rate of change, (c) a greater number of iterations per threshold value, and (d) a greater number of consecutive unsuccessful attempts prior to changing the threshold value. With respect to both forests (the Lincoln Tract and the Jones Tract), the results from employing threshold accepting (Tables 3 and 4) were generally what were expected. The highest quality solutions and the mean quality of solutions were located when using the more diligent search processes. Further, the lowest quality solutions were much worse when using the least diligent search processes than those located using the more diligent search processes. Interestingly, the variation in solution quality was relatively low when solving the two maximization problems, and relatively high when solving the two minimization problems, likely because the better solutions would ideally approach zero (thus the mean is lower) when solving the minimization problems.

Table 3.A comparison of the Lincoln Tract harvest scheduling outcomes from threshold accepting when using parameters that are reflective of the most and least diligent searches.

	Initial threshold value	Rate of threshold change	Successful iterations per threshold	Unsucc. Iterations per threshold	Highest quality solutiona	Mean quality solutiona	Lowest quality solutiona	Coefficient of variationb
1. Maximize net revenue / Unit restriction model of final harvest adjacency
Most diligent (3bβ)	5,000,000	0.9999	100	1000	43,448,190	43,019,276	42,240,028	0.51
Least diligent (1aα)	1,000,000	0.99	1	500	42,023,408	40,382,158	37,214,842	2.14
2. Maximize net revenue / Area restriction model of final harvest adjacency
Most diligent (3bβ)	5,000,000	0.9999	100	1000	43,831,050	43,460,759	43,123,240	0.34
Least diligent (1aα)	1,000,000	0.99	1	500	42,433,492	40,927,223	37,530,429	1.97
3. Minimize wood flow deviations / Unit restriction model of final harvest adjacency
Most diligent (3bβ)	5,000,000	0.9999	100	1000	3,225	127,118	5,792,699	460.13
Least diligent (1aα)	1,000,000	0.99	1	500	9,225	795,843	11,226,465	187.54
4. Minimize wood flow deviations / Area restriction model of final harvest adjacency
Most diligent (3bβ)	5,000,000	0.9999	100	1000	5,170	173,438	2,800,048	332.77
Least diligent (1aα)	1,000,000	0.99	1	500	15,544	978,801	11,249,614	186.91

^a Units are US dollars for the maximization problems and squared deviations of scheduled harvest volumes from a target harvest volume in each time period in the minimization problems

^b Percent, calculated as (standard deviation / mean) × 100

Table 4.A comparison of the Jones Tract harvest scheduling outcomes from threshold accepting when using parameters that are reflective of the most and least diligent searches.

	Initial threshold value	Rate of threshold change	Successful iterations per threshold	Unsucc. Iterations per threshold	Highest quality solutiona	Mean quality solutiona	Lowest quality solutiona	Coefficient of variationb
1. Maximize net revenue / Unit restriction model of final harvest adjacency
Most diligent (3bβ)	500,000	0.9999	100	1000	2,020,051	1,999,527	1,958,518	0.72
Least diligent (1aα)	200,000	0.99	1	500	2,011,828	1,975,977	1,924,821	0.91
2. Maximize net revenue / Area restriction model of final harvest adjacency
Most diligent (3bβ)	500,000	0.9999	100	1000	2,045,444	2,032,845	1,998,809	0.43
Least diligent (1aα)	200,000	0.99	1	500	2,022,865	1,958,843	1,795,534	2.61
3. Minimize wood flow deviations / Unit restriction model of final harvest adjacency
Most diligent (3bβ)	5,000,000	0.9999	100	1000	48	2,995	193,019	453.13
Least diligent (1aα)	1,000,000	0.99	1	500	1,130	104,503	1,982,576	263.70
4. Minimize wood flow deviations / Area restriction model of final harvest adjacency
Most diligent (3bβ)	5,000,000	0.9999	100	1000	37	842	2,469	58.19
Least diligent (1aα)	1,000,000	0.99	1	500	224	38,968	565,716	171.24

^a Units are US dollars for the maximization problems and squared deviations of scheduled harvest volumes from a target harvest volume in each time period in the minimization problems

^b Percent, calculated as (standard deviation / mean) × 100

When incorporating additional processes to the threshold accepting search process, the quality of the solutions generated generally improved. With the two maximization problems that were applied to Lincoln Tract, adding just 2-opt moves to the search process improved the best solution located even though the mean solution value declined in one of these cases (Table 5). Adding just the reversion process did not seem to improve the results in the two maximization problems, in fact the highest quality solution located, mean solution value, and lowest quality solution located were all lower than when reversion was not employed (simply using the most diligent search process parameters). This suggests that perhaps the search reversion rate employed was not the most appropriate for these problems. Further, with respect to the two maximization problems, when both 2-opt moves and search reversion are employed, an improvement in the highest quality solution was observed over simply using the most diligent search process parameters, yet the results were not as good as simply using 2-opt moves alone. With the two minimization problems, both 2-opt moves and search reversion, when employed separately, dramatically improved the results generated by the threshold search process. When employing both 2-opt moves and search reversion, the highest quality results were obtained.

Table 5.A comparison of the Lincoln Tract harvest scheduling outcomes from threshold accepting when using parameters that are reflective of the most diligent searches and further enhancements.

	Highest quality solutiona	Mean quality solutiona	Lowest quality solutiona	Coefficient of variationb
1. Maximize net revenue / Unit restriction model of final harvest adjacency
Most diligent (3bβ)	43,448,190	43,019,276	42,240,028	0.51
Most diligent (3bβ) + 2-opt movesc	43,552,305	43,288,773	42,723,856	0.31
Most diligent (3bβ) + search reversiond	43,359,417	42,370,812	40,802,585	1.02
Most diligent (3bβ) + 2-opt moves and search reversionc,d	43,507,808	42,819,229	41,754,987	0.71
2. Maximize net revenue / Area restriction model of final harvest adjacency
Most diligent (3bβ)	43,831,050	43,460,759	43,123,240	0.34
Most diligent (3bβ) + 2-opt movesc	43,882,583	42,924,888	37,816,493	3.49
Most diligent (3bβ) + search reversiond	43,705,603	42,992,724	41,676,410	0.85
Most diligent (3bβ) + 2-opt moves and search reversionc,d	43,782,704	43,131,960	40,889,511	1.05
3. Minimize wood flow deviations / Unit restriction model of final harvest adjacency
Most diligent (3bβ)	3,225	127,118	5,792,699	460.13
Most diligent (3bβ) + 2-opt movesc	155	2,607	21,650	83.11
Most diligent (3bβ) + search reversiond	744	110,949	4,036,578	454.98
Most diligent (3bβ) + 2-opt moves and search reversionc,d	24	926	5,970	98.59
4. Minimize wood flow deviations / Area restriction model of final harvest adjacency
Most diligent (3bβ)	5,170	173,438	2,800,048	332.77
Most diligent (3bβ) + 2-opt movesc	825	38,178	384,016	133.21
Most diligent (3bβ) + search reversiond	202	322,706	5,776,699	356.33
Most diligent (3bβ) + 2-opt moves and search reversionc,d	6	8,128	156,073	253.06

^a Units are US dollars for the maximization problems and squared deviations of scheduled harvest volumes from a target harvest volume in each time period in the minimization problems

^b Percent, calculated as (standard deviation / mean) × 100

^c For every fifty 1-opt moves, ten 2-opt moves are employed

^d For every fifty successful iterations, the search reverts back to the best solution stored in memory of the computer

When 2-opt moves and search reversion were applied to the maximization problems of the Jones Tract (Table 6), we observed that the addition of 2-opt moves improved the highest quality solutions located. In general, when either of these additional processes were added to the basic threshold accepting search process, the mean quality of solutions improved. Interestingly, the variation in solution quality was much lower when search reversion was employed compared to simply using the most diligent search process parameters or employing 2-opt moves alone. In the two minimization problems, both 2-opt moves and search reversion, when employed separately, dramatically improved the results generated by the threshold search process. And again, when employing both 2-opt moves and search reversion, the highest quality results were obtained. In fact, the optimal solutions to the two minimization problems were located using the threshold accepting search process supplemented with both 2-opt moves and search reversion.

Table 6.A comparison of the Jones Tract harvest scheduling outcomes from threshold accepting when using parameters that are reflective of the most diligent searches and further enhancements.

	Highest quality solutiona	Mean quality solutiona	Lowest quality solutiona	Coefficient of variationb
1. Maximize net revenue / Unit restriction model of final harvest adjacency
Most diligent (3bβ)	2,020,051	1,999,527	1,958,518	0.72
Most diligent (3bβ) + 2-opt movesc	2,020,051	2,005,953	1,967,365	0.59
Most diligent (3bβ) + search reversiond	2,019,279	1,976,534	1,898,641	1.01
Most diligent (3bβ) + 2-opt moves and search reversion	2,019,460	1,984,325	1,924,852	1.06
2. Maximize net revenue / Area restriction model of final harvest adjacency
Most diligent (3bβ)	2,045,444	2,032,845	1,998,809	0.43
Most diligent (3bβ) + 2-opt movesc	2,045,594	2,022,447	1,857,429	2.06
Most diligent (3bβ) + search reversiond	2,043,983	2,024,873	1,988,226	0.53
Most diligent (3bβ) + 2-opt moves and search reversion	2,045,581	2,028,689	1,953,502	0.61
3. Minimize wood flow deviations / Unit restriction model of final harvest adjacency
Most diligent (3bβ)	48	2,995	193,019	453.13
Most diligent (3bβ) + 2-opt movesc	4	81	373	76.30
Most diligent (3bβ) + search reversiond	10	329	4,253	134.49
Most diligent (3bβ) + 2-opt moves and search reversion	0e	12	93	119.56
4. Minimize wood flow deviations / Area restriction model of final harvest adjacency
Most diligent (3bβ)	37	842	2,469	58.19
Most diligent (3bβ) + 2-opt movesc	7	1,128	15,792	183.27
Most diligent (3bβ) + search reversiond	2	121	1,208	138.96
Most diligent (3bβ) + 2-opt moves and search reversion	0f	46	1,467	299.14

^a Units are US dollars for the maximization problems and squared deviations of scheduled harvest volumes from a target harvest volume in each time period in the minimization problems

^b Percent, calculated as (standard deviation / mean) × 100

^c For every fifty 1-opt moves, ten 2-opt moves are employed

^d For every fifty successful iterations, the search reverts back to the best solution stored in memory of the computer

^e 0.05

^f 0.02

The mixed integer programming objective function values to the maximization problems applied to the Lincoln Tract were $43,581,416.80 (URM) and $43,882,583.90 (ARM). The threshold accepting search process located the ARM optimal solution twice (in 1% of the solutions generated) using 3bβ parameters and 2-opt moves employed. The very best threshold accepting solution for the URM problem was $56,161.80 (0.13%) from the URM optimal solution using 3bα● parameters without the use of 2-opt moves or search reversion. The mixed integer programming objective function values to the maximization problems applied to the Jones Tract were $2,020,050.97 (URM) and $2,045,594.00 (ARM). The threshold accepting search process located the URM optimal solution once using 3bβ parameters without the enhancement of 2-opt moves or search reversion, and once using 3bβ parameters with the enhancement of 2-opt moves. The threshold accepting search process located the ARM optimal solution one time using using 3bβ parameters with the enhancement of 2-opt moves. In total, for the maximization problems applied to both case study forests, threshold accepting located the optimal solution 5 times out of 31,200 attempts using the sets of search parameters studied here (0.016% of the solutions generated). However, 47.1% of the solutions generated by threshold accepting were within 1% of the global optimum values for the maximization problems applied to both case study forests.

The mixed integer programming objective function values to the minimization problems applied to the Lincoln Tract were 0.90 (URM) and 0.45 (ARM). The very best threshold accepting solution for the URM problem was 24.0 using 3bβ parameters and enhanced with 2-opt moves and search reversion. Similarly, the very best threshold accepting solution for the ARM problem was 5.8 using 3bβ parameters and enhanced with 2-opt moves and search reversion. The mixed integer programming objective function values to the minimization problems applied to the Jones Tract were 0.05 (URM) and 0.02 (ARM). The very best threshold accepting solution for the URM problem was 0.32 using 3bβ parameters and enhanced with 2-opt moves and search reversion. The threshold accepting search process located the ARM optimal solution one time using using 3bβ parameters with the enhancement of 2-opt moves and search reversion. In total, for the minimization problems applied to both case study forests, threshold accepting located the optimal solution 1 time out of 15,600 attempts using the sets of search parameters studied here.

Tests of Significant Differences

Most paired tests of significant differences (p < 0.05) amongst sets of outcomes (200 samples per parameter set) indicated that the objective function values of the sets were statistically significantly different, with a few exceptions. From the vast array of paired comparisons amongst two different sets of solutions generated by the threshold accepting heuristic, some general trends were evident.

• 1. Within a forest management problem, when the heuristic search was applied to either of the two case study forests using (a) the same rate of change in the threshold, (b) the same number of iterations per threshold, and (c) the same number of unsuccessful iterations per threshold, yet the search was initiated from different threshold levels, there were generally no statistically significant differences (p = 0.05) in the quality of the objective function values within sets of solutions generated.

• 2. Within a forest management problem, when the heuristic search was applied to either of the two case study forests using (a) the same initial threshold level, (b) the same rate of change in the threshold, (c) only 1 iteration per threshold, yet a different number of unsuccessful iterations per threshold were assumed (500 or 1000), there were generally no statistically significant differences (p = 0.05) in the quality of the objective function values within sets of solutions generated. This was not the case when 100 iterations per threshold were assumed.

• 3. Within a forest management problem, when the heuristic search was applied to either of the two case study forests using (a) the same initial threshold, (b) the same number of iterations per threshold, (c) the same number of unsuccessful iterations per threshold, yet a different rate of change in the threshold, there were generally statistically significant differences (p = 0.05) in the quality of the objective function values within sets of solutions generated.

• 4. Within a forest management problem, when the heuristic search was applied to either of the two case study forests using (a) the same initial threshold, (b) the same rate of change in the threshold value, (c) the same number of unsuccessful iterations per threshold, yet different iterations per threshold (1 vs. 100), there were generally statistically significant differences (p = 0.05) in the quality of the objective function values within sets of solutions generated.

• 5. Within the maximization problems, when the heuristic search was applied to either of the two case study forests using (a) a moderate rate of change (0.999), (b) different initial threshold levels, (c) different number of unsuccessful iterations per threshold, and (c) where the iterations per threshold was assumed to be 1000, there were generally no statistically significant differences (p = 0.05) in the quality of the objective function values when these alternatives were compared with the most diligent searches using 3bβ parameters enhanced with 2-opt moves and search reversion. However, generally, statistically significant differences were observed when these sets of solutions were compared to searches using 3bβ parameters enhanced with 2-opt moves alone. Furthermore, these findings were not observed in the minimization problems.

• 6. Within the minimization problems, in nearly all instances of different search parameters, the sets of solutions from the three alternatives that included enhancements (2-opt moves, search reversion) were statistically significant different (p = 0.05) with respect to the quality of the objective function values within other 36 sets of solutions generated.

DISCUSSION

When applied to the forest management problems described in this paper, threshold accepting acted as a relatively fast search heuristic. In general, with each run of the heuristic search process there is a hill-climbing phase, an adjustment phase, and a fine-tuning or steady-state phase (Fig. 10) (Bettinger et al., 1997). However, specific search parameters can affect the amount of time needed to complete the search process. For example, there seemed to be a 10-fold increase in the amount of time required as the assumed rate of change increased in magnitude. When the rate of change of the threshold value was 0.99, the general amount of time required to complete the search process was about 6 seconds per run. As the rate of change increased to 0.999, about one minute per run was required, and as the rate of change increased to 0.9999 about 10 minutes per run were required. It is acknowledged that these observations of performance can improve had the search process been coded in another computer language or enabled on a computer with faster CPU speed. The speed at which exact methods solved the maximization problems suggests that for forest harvest scheduling cases where the adjacency constraints can easily be developed and the wood flow constraints are not too narrow, this may be the more appealing option. When there may be difficulty in developing adjacency constraints (particularly in ARM problems with multiple time period green-up assumptions) or when there is a desire for very narrow wood flow deviations from period to period, the heuristic search process may be the more appealing option. The reason, for example, that the minimization problems required an extensive amount of time when branch and bound methods were employed was the desire for a harvest schedule that contained no deviation in the wood flow from period to period.

Fig. 10.

General behavior of the threshold accepting heuristic search process when it is applied to the Lincoln Tract and forest management problem 1 (maximization of discounted net revenue and unit restriction adjacency of final harvests).

Although a large number of heuristic search solutions were developed for this research effort, we could have explored different assumptions regarding the implementation of 2-opt moves and we could have explored different assumptions regarding the search reversion rate (as has been shown in Bettinger et al. (2015)). Should these activities eventually be pursued, the most appropriate combination of 1-opt and 2-opt moves might be assessed, along with the most appropriate search reversion rate. We hypothesize that the most appropriate of these options, in combination with the normal heuristic search parameters, may be problem-specific. Since the intent of this research effort was to illustrate the general quality of solutions generated by the threshold accepting heuristic search process as parameter values are altered, a deeper analysis of enhancements is left open for future research efforts.

While the optimal solution to several of the problems investigated here was located, Winkler and Fang (1997) noted that beyond some acceptable amount of computing effort, further increases in the number of iterations of local search algorithms generally will not improve the overall quality of solutions produced. And so, the hope of locating higher quality solutions may be just that (wishful thinking) when using parameter sets that are unable to allow a diligent search through a solution space. Enhancements such as those described here (2-opt moves and search reversion) that are embedded in s-metaheuristic search processes may overcome some of these difficulties (Bettinger and Boston, 2017). Further, slowing the rate of change of the threshold value or increasing the number of consecutive unsuccessful moves needed before changing the threshold value both act to increase the number of iterations of this search process. To be consistent with the thoughts of Winkler and Fang (1997) all sets of solutions less than the most diligent search would seem to not employ an acceptable amount of computing effort as the most diligent searches led to improved objective function values.

There are other methods for adjusting the manner in which threshold accepting conducts a search. For example, Lin et al. (1995) described how the search might be limited to only promising areas of the solution space, rather than consideration of all random moves from the current solution to the proposed solution. The reversion process that is described in this work is similar and serves to also limit the search to promising areas of the solution space. Further, Lin et al. (1995) described methods for adapting the threshold value based on recent search performance, allowing the threshold to not only decrease as the search progresses, but also to increase again when many local minima and maxima are recognized during a short period of search time. This form of adjustment to the acceptance criteria is similar in concept to what might be employed in a demon algorithm (Creutz, 1983; Wood and Downs, 1998), which is a s-metaheuristic that utilizes an acceptance rule that varies in magnitude depending on recent search history. Other means for adjusting a search process have been investigated, for example, Nissen and Paul (1995) suggested increasing the number of iterations allowed as the threshold declines, since the ability to locate improved solutions and escape from local optima may be more difficult when the threshold value decreases. Further, Frausto-Solis et al. (2021) suggested a form of a ruin and recreate process (as they note, a chaos perturbation) to re-start a threshold accepting search process from a different place in the solution space. Even further, Dhouib et al. (2010) suggested that to diversity a search, one might add memory to the search process, through the use of a tabu state assigned only to moves that lead to a lower quality objective function value. These opportunities to perhaps improve the manner in which a threshold accepting search process addresses a forest harvest scheduling problem are left to others to pursue.

As a final statement on the usefulness of a threshold accepting in forest harvest scheduling efforts, both Zhu and Bettinger (2008) and Li et al. (2010) have suggested that this form of s-metaheuristic search process could be combined with a different form of s-metaheuristic search process (tabu search) to improve the quality of forest harvest schedules that have wood flow and harvest adjacency constraints. In these cases, a threshold accepting search process may be beneficial for the initial hill-climbing phase of the search due to the speed at which the process moves an initial randomly developed (and poor quality) forest harvest schedule to a moderately well-valued schedule. A tabu search heuristic could then the employed upon the best solution provided by threshold accepting to adjust and refine the quality of the harvest schedule. As tabu search needs to assess a number of potential neighboring solutions during each iteration of the model prior to selecting one, it is effectively much slower in operation than threshold accepting. Yet the value of adding tabu search may be in the deterministic nature that it employs to adjust and refine the threshold accepting solution. As the potential value of the combined metaheuristic approach has been described in both Zhu and Bettinger (2008) and Li et al. (2010), and since this research effort focuses simply on the behavior of threshold accepting, readers are encouraged to visit the other works to understand more fully the value of the combined approach.

CONCLUSIONS

From the analysis provided in this research effort it is evident that high quality solutions to forest harvest scheduling problems can be produced as a threshold accepting heuristic search process more diligently navigates through the solution space. This diligence involves allowing higher amounts of unsuccessful attempts allowed before the threshold value changes, a slower rate of change in the threshold value, the use of 2-opt moves to navigate through the neighborhood, and periodic reversion to the best solution held in memory to intensify the search in high quality areas of the solution space. Significant differences in the quality of the forest harvest schedules developed using threshold accepting, as reflected by the objective function values, were observed when the number of iterations per threshold varied and when the rate of change in the threshold value varied. Although the number of unsuccessful iterations per threshold and the initial threshold level are important assumptions for the heuristic search process, in certain cases, no significant differences in solution quality were observed when these varied. As a general rule, particularly when wood flow deviations are narrow (as in the minimization problem investigated here), enhancements to threshold accepting (2-opt moves and search reversion) resulted in higher quality forest harvest scheduling solutions when used in conjunction with the search parameters that represented the most diligent search.

LITERATURE CITED

• Akbulut, R., Bettinger, P., Ucar, Z., Obata, S., Boston, K. and Siry, J. (2017) Spatial forest plan development using heuristic processes seeded with a relaxed linear programming solution. Forest Sci. 63 : 518–528
• Albus, J.S. (1981) Brains, behavior, and robotics. BYTE Publications, Inc., Peterborough, New Hampshire, 352 pp
• Bettinger, P. and Boston, K. (2008) Habitat and commodity production trade-offs in coastal Oregon. Socio-Econ. Plan. Sci. 42 : 112–128
• Bettinger, P. and Boston, K. (2017) Forest planning heuristics - Current recommendations and research opportunities for s-metaheuristics. Forests. 8 : Article 476
• Bettinger, P., Boston, K. and Sessions, J. (1999) Intensifying a heuristic forest harvest scheduling search procedure with 2-opt decision choices. Can. J. Forest Res. 29 : 1784–1792
• Bettinger, P., Boston, K., Siry, J.P. and Grebner, D.L. (2017) Forest management and planning, second edition. Academic Press, London, 349 pp
• Bettinger, P., Demirci, M. and Boston, K. (2015) Search reversion within s-metaheuristics: Impacts illustrated with a forest planning problem. Silva Fenn. 49 : article id 1232
• Bettinger, P., Graetz, D., Boston, K., Sessions, J. and Chung, W. (2002) Eight heuristic planning techniques applied to three increasingly difficult wildlife planning problems. Silva Fenn. 36 : 561–584
• Bettinger, P., Johnson, D.L. and Johnson K.N. (2003) Spatial forest plan development with ecological and economic goals. Ecol. Model. 169 : 215–236
• Bettinger, P., Sessions, J. and Boston, K. (1997) Using Tabu search to schedule timber harvests subject to spatial wildlife goals for big game. Ecol. Model. 94 : 111–123
• Brodie, J.D. and Sessions, J. (1991) The evolution of analytic approaches to spatial harvest scheduling. In: Buford, M.A. (comp.) Proceedings of the 1991 Symposium on Systems Analysis in Forest Resources, General Technical Report SE-74. U.S. Department of Agriculture, Forest Service, Southeastern Forest Experiment Station, Asheville, North Carolina: 187–191
• Caro, F., Constantino, M., Martins, I. and Weintraub, A. (2003) A 2-opt tabu search procedure for the multiperiod forest harvesting problem with adjacency, greenup, old growth, and even flow constraints. For. Sci. 49 : 738–751
• Coulter, E.D., Sessions, J. and Wing, M.G. (2006) Scheduling forest road maintenance using the analytic hierarchy process and heuristics. Silva Fenn. 40 : 143–160
• Creutz, M. (1983) Microcanonical Monte Carlo simulation. Phys. Rev. Lett. 50 : 1411–1414
• Dhouib, S., Kharrat, A. and Chabchoub, H. (2010) A multi-start threshold accepting algorithm for multiple objective continuous optimization problems. Int. J. Numer. Meth. Eng. 83 : 1498–1517
• Dueck, G. and Scheuer, T. (1990) Threshold accepting: A general purpose optimization algorithm appearing superior to simulated annealing. J. Comput. Phys. 90 : 161–175
• Frausto-Solis, J., Hernández-Ramírez, L., Castilla-Valdez, G., González-Barbosa, J.J. and Sánchez-Hernández, J.P. (2021) Chaotic multi-objective simulated annealing and threshold accepting for job shop scheduling problem. Math. Comput. Appl. 26 : Article 8
• Golden, B.L. and Alt, F.B. (1979) Interval estimation of a global optimum for large combinatorial problems. Nav. Res. Logist. Q. 26 : 69–77
• Hanzlik, E.J. (1922) Determination of the annual cut on a sustained basis for virgin American forests. J. For. 20 : 611–625
• Killen, J. (2021) Mathematical programming methods for geographers and planners. Routledge, New York, 386 pp
• Li, R., Bettinger, P. and Boston, K. (2010) Informed development of meta heuristics for spatial forest planning problems. The Open Oper. Res. J. 4 : 1–11
• Lin, C.K.Y., Haley, K.B. and Sparks, C. (1995) A comparative study of both standard and adaptive versions of threshold accepting and simulated annealing algorithms in three scheduling problems. Eur. J. Oper. Res. 83 : 330–346
• LINDO Systems Inc. (2024) Extended Lingo/Win64 release 21.0.32. LINDO Systems Inc., Chicago, Illinois
• Los, M. and Lardinois, C. (1982) Combinatorial programming, statistical optimization and the optimal transportation network problem. Transp. Res. - Part B. 16B : 89–124
• Marshall, H., Boston, K. and Sessions, J. (2003) Harvest scheduling with aggregation adjacent constraints: Threshold acceptance approach. In: Briggs, D. (chair) Precision forestry, proceedings of the second international precision forestry symposium. University of Washington, College of Forest Resources, Institute of Forest Resources, Seattle: 131–135
• Meyer, H.A. (1952) Structure, growth, and drain in balanced even-aged forests. J. Forest. 50: 85–95
• Murray, A.T. (1999) Spatial restrictions in harvest scheduling. For. Sci. 45 : 45–52
• Nissen, V. and Paul, H. (1995) A modification of threshold accepting and its application to the quadratic assignment problem. OR Spektrum. 17 : 205–210
• Perea, C., Alcala, J., Yepes, V., Gonzalez-Vidosa, F. and Hospitaler, A. (2008) Design of concrete bridge frames by heuristic optimization. Adv. Eng. Softw. 39 : 676–688
• Pukkala, T. and Heinonen, T. (2006) Optimizing heuristic search in forest planning. Nonlinear Anal.- Real. 7 : 1284–1297
• Romanycia, M.H. and Pelletier, F.J. (1985) What is a heuristic? Comput. Intell. 1 : 47–58
• Siedentopf, J. (1995) An efficient scheduling algorithm based on threshold accepting. Institute of Production Management and Industrial Information Management, University of Leipzig, Leipzig, Germany. Report No. 16
• Talbi, E.-G. (2009) Metaheuristics from design to implementation. John Wiley & Sons, Inc., Hoboken, New Jersey, 593 pp
• Tarantilis, C.D. and Kiranoudis, C.T. (2002) A threshold-accepting metaheuristic method for scheduling the operations of dehydration plants. Dry. Technol. 20 : 1143–1160
• Tarantilis, C.D., Kiranoudis, C.T. and Vassiliadis, V.S. (2004) A threshold accepting metaheuristic for the heterogeneous fixed fleet vehicle routing problem. Eur. J. Oper. Res. 152 : 148–158
• Winkler, P. and Fang, K.-T. (1997) Application of threshold-accepting to the evaluation of the discrepancy of a set of points. SIAM J. Numer. Anal. 34 : 2028–2042
• Wood, I. and Downs, T. (1998) Demon algorithms and their application to optimization problems. In: Simpson, P.K. (general chair) IEEE world congress on computational intelligence. Institute of Electrical and Electronics Engineers, Inc., Piscataway, New Jersey: 1661–1666
• Zhu, J. and Bettinger, P. (2008) Estimating the effects of adjacency and green-up constraints on landowners of different sizes and spatial arrangements located in the southeastern U.S. Forest Policy Econ. 10 : 295–302