Train corridor scheduling process including various cost functions associated with railway operations

[0002] This invention relates to a process for scheduling the movement of trains over a rail corridor having a plurality of sidings or parallel tracks with crossover switches.

[0003] A rail corridor is a collection of tracks and sidings connecting two rail terminal areas. An example of a rail corridor 8 is shown in FIG. 1, showing a single main track 10 and three sidings 20. The western end of the rail corridor is on the left side of FIG. 1 and the eastern end on the right.

[0004] Scheduling rail transportation on a rail corridor is particularly complex as compared to highway, water, or air transportation. Trains using a single track traveling in opposite directions (i.e., a meet) or trains traveling in the same direction (i.e., a pass) must meet in the vicinity of a siding so that one train can be sided to let the other pass. Alternatively, if there exists a double main line with crossover switches, one train can be switched to the second main line to allow the other train to pass. Also, when such meets or passes occur at a siding, the siding chosen must be long enough to accommodate the train to be sided, and the train to be sided must arrive at the siding and have sufficient time to pull onto the siding before the passing train arrives at the siding.

[0005] The railroad must earn revenue from its transportation operations, and some of this revenue is generally at risk if trains cannot deliver freight on time. The destination time of the trains must be managed insofar as possible to prevent late penalties incurred by the railroad. Therefore scheduling trains across a rail corridor involves arranging meets and passes as required for all trains, and while also meeting the schedule for each train so that they all arrive, on time, at the end of the corridor.

[0006] Commercially applied scheduling processes attempted to date have been based on paradigms which involve simulation with branch and bound techniques to find a conflict-free schedule. Since a branch and bound process must sort through many binary choices as it proceeds toward a solution, these techniques are slow, and do not take advantage of quantitative relationships that can be adduced from the scheduling context.

[0007] Additionally, the prior art technique search processes actually become more complex and take longer to arrive at a solution as the number of sidings in the rail corridor increases. This is due to the search algorithms that form the basis for these prior art techniques. More sidings requires the search algorithm to search through and consider more choices before arriving at an optimum solution. As will be shown below, the technique of the present invention overcomes this disadvantage. Since the present invention calculates a cost function where each siding represents a lower cost, having more sidings will make it easier for the algorithm to identify the optimal (i.e. minimal) cost.

[0008] One prior art technique uses quantitative information such as train speed, destination, and time of departure as discrete variables in an artificial intelligence based system. The artificial intelligence process involves rules that are used to search through the trial cases until the best case is found. In addition to the considerable time taken by an artificial intelligence system to optimize a solution, it is also known that a slight change to the initial conditions may produce a significantly different result. In any case, a slight change to the initial conditions will require a new and lengthy computation to find the optimum solution. A commercial product referred to as The Movement Planner, offered by G E-Harris Railway Electronics L.L.C. of Melbourne, Fla., implements such an artificial intelligence solution.

[0009] As can be seen, the total set of parameters for scheduling a corridor can be large, and of both discrete and continuous types. Generally, a cost function based on these parameters can be formulated, and then some method of search is executed that will reduce the cost and/or find a feasible schedule for the subject trains. But, the presence of discrete variables in the search space prevents or greatly complicates the application of any “hill-climbing” search processes based on the use of gradients

[0010] Cost functions that are everywhere differentiable have the advantage over prior art artificial intelligence solutions of being amenable to gradient-based minimization algorithms that do not have to accommodate the difficulties that arise in discrete or partially discrete search spaces. The present invention is a process whereby a rail corridor and the train schedule along that corridor can be characterized by a differentiable (i.e., continuous) cost function, so that a search process based on differentiation may be applied to scheduling train activity in the corridor.

[0011] The present invention is an analytical process for scheduling trains across a corridor that is driven by a cost function to be minimized, where the cost function is a continuous and differentiable function of the scheduling variables. The present invention is an improvement over the prior art cost functions that include discrete variables and thus are not differentiable everywhere. The present invention will permit the use of search processes relying on gradients, and as such, will converge to solutions much more quickly than the prior art scheduling processes involving simulation, or searching through discrete options.

[0012] The corridor scheduling process of the present invention involves three steps for identification of the optimum schedule. After an acceptable differentiable cost function is derived, the first step is the gradient search process wherein the gradient of the differentiable cost function is determined. The cost function is a sum of individual localizer functions. For each pair of trains in the corridor that might intersect, using the localizer function, the intersection point is identified as having a high value if the train trajectories do not intersect near a siding and lower values as the intersection point moves toward any siding. The gradient process may not move all intersection points precisely to the center of sidings dependent upon the selected threshold value and parametric values of the localizer function. Instead, the gradient process varies train departure times so that the set of all intersection points of trains are moved nearer to sidings. The second phase of the process simply moves the points precisely to the centers of sidings, selects which train to side, and computes exact arrival and departure times for the trains at the siding to assure the physical integrity of the meet. In order to center the intersection points at sidings and side specific trains, the speeds of the individual trains must be modified. This is accomplished during the second step of the scheduling process.

[0013] The third step maintains the proper siding relationships between any two meeting trains, as determined in step two, but allows the meet time to vary in an effort to assure that no train exceeds an upper speed limit. This final phase is again a gradient search process applied to all of the meet points determined in the second step.

[0014] The present invention can be more easily understood, and the further advantages and uses thereof more readily apparent, when considered in view of the description of the preferred embodies and the following figures. Identical reference characters in the figures identical components of the invention.

[0015] FIG. 1 illustrates of a simple rail corridor;

[0016] FIG. 2 is a string diagram illustrating the corridor scheduling problem in terms of intersecting lines;

[0017] FIG. 3 is a flow chart for the corridor scheduling process of the present invention;

[0018] FIG. 4 illustrates the basic geometry of train trajectories;

[0019] FIG. 5 is a graph of the basic sigmoid function;

[0020] FIG. 6 illustrates the use of sigmoid sums to discriminate an interval;

[0021] FIG. 7 illustrates the construction of a localizer function from sigmoid functions;

[0022] FIG. 8 illustrates an example of a localizer function for two sidings;

[0023] FIGS. 9A and 9B show the modification of a localizer function to account for corridor endpoints;

[0024] FIG. 10 illustrates the necessary geometry to achieve a balanced localizer function;

[0025] FIGS. 11A, 11B, and 11C illustrate a technique for approximating the economic penalty function;

[0026] FIG. 12 shows a penalty term function for early departure of a train;

[0027] FIG. 13 is an initial infeasible string graph schedule for twelve trains;

[0028] FIG. 14 is a string graph for trains of FIG. 13 after a gradient search of the present invention;

[0029] FIG. 15 shows the process whereby intersection points are moved to a siding center;

[0030] FIG. 16 shows moving the first intersection point to a siding center;

[0031] FIG. 17 illustrates the process of speed adjustments to center all meets;

[0032] FIGS. 18A and 18B through FIGS. 24A and 24B illustrate certain infeasibilities created by centering meets on sidings and the resolution thereof;

[0033] FIGS. 19A and 19B illustrates the two types of siding conflict;

[0034] FIGS. 20A and 20B illustrate the resolution of certain siding conflicts;

[0035] FIGS. 21A and 21B illustrate the “unresolvable” siding conflict;

[0036] FIGS. 22A through 22D illustrate resolution of both types of siding conflicts;

[0037] FIGS. 23A through 23E show the cases for downward resolvable siding conflicts;

[0038] FIGS. 24A and 24B show the resolution of upward-resolvable siding conflicts;

[0039] FIG. 25 illustrates train trajectories represented as broken line segments;

[0040] FIG. 26 is an evaluation of the train trajectory vector;

[0041] FIG. 27 shows an adjustment of train trajectory to accommodate siding delays;

[0042] FIG. 28 shows siding details for a westbound sided train;

[0043] FIG. 29 illustrates siding details for eastbound passing trains;

[0044] FIG. 30 is a complete string graph adjusted for centered meets and train sidings; and

[0045] FIGS. 31 and 32 are flow charts illustrating algorithms implemented by the present invention.

[0046] Before describing in detail the train corridor scheduling process in accordance with the present invention, it should be observed that the present invention resides primarily in a novel scheduling process algorithm and not in the particular detailed configurations thereof.

[0047] The traditional method of graphic depiction of a train schedule for a rail corridor is referred to as a string graph as shown in FIG. 2. This string graph represents a time-distance graph of train movement in the corridor depicted in FIG. 1. The horizontal axis represents time, (i.e., a fixed window of time) and the vertical axis represents distance, with the point at the origin of the graph being the western end of the corridor, and the point at the top being the eastern end of the corridor. The width of the graph represents the period of interest in which the trains will be scheduled. Lines on the graph sloping one way represent traffic in one direction across the corridor, while lines slopping in the opposite direction represent oppositely-directed traffic. Only the position of the engine is shown. The horizontal bars across the graph, bearing reference to character 20, correspond to the siding locations.

[0048] The invention as presented herein is described in conjunction with a single-rail corridor with sidings. But those skilled in the art will recognize that it can be easily extended to multiple track main lines with cross-over switches between the main lines.

[0049] The essential criterion for an acceptable schedule, as expressed in terms of the string graph of FIG. 2, is that any two train trajectories (lines) on the graph must intersect at a siding 20. If their meets are at sidings, then in addition, a choice has to be made as to which train to side.

[0050] Note that, unless all of the intersecting lines actually intersect within the sidings 20, the schedule is infeasible. Assuming, for the nonce, that all train speeds will be fixed, the departure times for the trains can be adjusted in order to move the train lines about and attempt to place all intersection points over the sidings 20. In another embodiment of the present invention, it would be possible, as well, to vary train speeds, which would change the slopes of the train trajectory lines, in order to place intersection points over the sidings 20. In yet another embodiment, both speeds and departure times can be varied simultaneously to find a feasible meet/pass plan for the trains.

[0051] The process to be described herein treats the corridor scheduling problem as a geometry problem, rather than directly as a scheduling problem, as suggested by the prior art. It does so by providing a mechanism by which train trajectory lines are moved under control of a gradient-search process based on a differentiable cost function in a manner that moves the intersection points to or close to established sidings.

[0052] The search process of the present invention permits variation of speeds and departure times, separately or jointly, and will use an everywhere differentiable cost function that takes on lower values as the schedule approaches feasibility. Because the cost function is everywhere differentiable, an iterative, gradient-search method can be applied that assures that the successive schedules found by the search process in fact converge to a conflict-free result.

[0053] Moreover, it is possible to include, in another embodiment of the present invention, the constraint that a siding must be longer than a train to be sided on it. It is further possible to include, in yet another embodiment, the economic costs incurred by adjusting train schedules. In other embodiments, constraints on maximum train speed and the early departure of trains can also be considered.

[0054] It will be appreciated by those skilled in the art that although FIG. 2 illustrates a situation with three sidings and three trains traveling in each direction, the technique of the present invention can be easily extended to any number of trains operating in each direction and any number of sidings on the rail corridor. The concepts of the present invention can also be extended to a rail corridor with more than one main line and crossover switches between the main line tracks. The present invention can be applied to any rail corridor where one train can be switched to another track when a meet or pass with another train occurs.

[0055] The scheduling of trains must first be feasible, but in addition, there may be choices as to which trains to side or the order to run trains, which helps to assure that economic penalties will not be incurred or, failing that, will at least be ameliorated.

[0056] Process 30 for obtaining both schedule feasibility and economic acceptability may consist of a number of steps as shown in FIG. 3. First, at step 31, an initial prearrangement of the trains is done, establishing their order of entry into the corridor. At this point, the train order is based solely on due times, (represented as an input to step 31 from block 32) with no analysis as to the corridor capacity or specific departure times. At step 33, an initial schedule for the trains is determined; there are several numerical optimization techniques that may be applied here. See for example, Numerical Optimizationby Jorge Nacedad and Stephen J. Wright; Springer, New York 1999; ISBN 0-387-98793-2.

[0057] This initial schedule is input to the gradient search process, step 34, to be discussed below, which minimizes schedule infeasibility. In another embodiment the gradient search process can also minimize economic penalties incurred by the railroad for the late arrival of trains and give due consideration to maximum train speeds, early departure times and siding lengths. The gradient search adjusts train departure times (i.e., the time the train enters the corridor) and/or speeds so that meets occur near sidings. The process 30 loops through siding choice step 38 and the conflicts decision step 36 until all train intersections are placed at or near sidings on the rail corridor by adjusting the speed and/or departure time (i.e., the time the train enters the corridor) of the trains traversing the corridor.

[0058] The decisions made at step 38 as to which train to side for each pair of trains meeting at a siding may be driven by considerations of relative economic cost due to the delays created by siding one train versus another train. This siding decision process represents another embodiment of the present invention and will be discussed further below.

[0059] Once the siding decisions are made, some of the trajectories (those for sided trains) on the string graph (FIG. 2) will become broken lines, (representing infeasible meets) which may cause new schedule infeasibilities for some train trajectories. At this point, the gradient search can again be applied, but only to the subset of subtrajectories that have been driven into infeasible meets. Multiple passes through the gradient search step 34 and siding decision process step 38 should bring the schedule to complete feasibility.

[0060] FIG. 4 characterizes the train trajectories as lines based on the initial departure times (moment of entry into the corridor) and the train speed. In FIG. 4, the bottom of the vertical axis represents the west end of the corridor, and the positive direction along that axis corresponds to eastbound travel. The time window of interest for travel in the corridor begins at time d₀, and the length of the corridor is denoted by L.

[0061] FIG. 4 focuses on characterizing one eastbound train and one westbound train, respectively T_iand T_j, with corresponding trajectories labeled L_iand L_j. s_i, s_jdenote the speeds and d_i, d_jdenote the departure times of the trains T_iand T_jrespectively. The departure time of a train is the time at which it enters the corridor: for an eastbound train, that corresponds to a point situated on the horizontal axis of FIG. 4, (i.e., t=0) and for a westbound train, that corresponds to a point located on the horizontal line y=L.

[0062] Then for train trajectory L_i(eastbound), we may express the relationship between coordinates for any point on the line in the form $\begin{array}{cc}\frac{y-0}{t-d_i}=s_i,\mathrm{or}y=s_{i}t-s_i{d}_i.& \text{(3-1)}\end{array}$

[0063] For train T_j(westbound), the form of trajectory L_jcan be likewise expressed as $\begin{array}{cc}\frac{y-L}{t-d_j}=-s_j,\mathrm{or}y=-s_{j}t+s_j{d}_j+L.& \text{(3-2)}\end{array}$

[0064] We can write equations of identical form for both eastbound and westbound trains by writing

[0065] y=s_it−s_id_i+θ_iL  (3-3),

[0066] where the speed of westbound trains by convention will be the negative of the train's actual speed, and $\begin{array}{cc}{\theta }_i=\{\begin{array}{c}0\mathrm{if}\mathrm{train}{T}_i\mathrm{is}\mathrm{eastbound}\\ 1\mathrm{if}\mathrm{train}{T}_i\mathrm{is}\mathrm{westbound}.\end{array}& \text{(3-4)}\end{array}$

[0067] This form of a linear equation (3-3) is not the usual form directly in terms of slope and intercept, but in this analysis train speeds and departure times will be varied and the form of Equation 3-3 has the advantage of expressing the train trajectories explicitly in terms of speeds and departure times.

[0068] The objective of the present invention is to determine the coordinates of intersection points (t_ij, y_ij) for pairs of train trajectories, and move these intersection points to sidings. For trains T_iand T_j, the solution for the trajectory intersection point is (t_ij, y_ij), where $\begin{array}{cc}{t}_{\mathrm{ij}}=\frac{s_i{d}_i-s_j{d}_j+\left({\theta }_j-{\theta }_i\right)L}{s_i-s_j},\mathrm{and}& \text{(3-5)}\\ y_{\mathrm{ij}}=\frac{s_i{s}_j\left(d_i-d_j\right)+\left(s_i{\theta }_j-s_j{\theta }_i\right)L}{s_i-s_j}.& \text{(3-6)}\end{array}$

[0069] (t_ij, y_ij) is derived by equating equation (3-1) and (3-2) (after making the notation change suggested by equation (3-3)).

[0070] This characterization of the intersection point applies to intersections of like-directed or oppositely-directed trains, so that the analysis to be developed concerning intersection points will adjust train trajectories involving both meets and passes.

[0071] To this point the train scheduling problem has been abstracted to a context of moving intersecting lines about until all intersection points are within certain ranges (the siding bars 20 in FIG. 2). When all intersection points that are within the rectangle (representing the corridor 8) are also within the sidings 20, we have obtained a feasible schedule.

[0072] It is an objective to obtain a feasible schedule using a search process that minimizes a cost function, and in the preferred embodiment, the preferred cost function will have a high value if any intersection point is outside a siding bar, and a low value if and only if all intersection points are within the siding bars. Intersection points entirely outside the graph are not considered; the corridor and scheduling period are considered to be co-extensive with the graph.

[0073] Let a function of a single y_ijwith this cost function property be a localizer function, and construct such a localizer function using the sigrnoid function as a basis. The preferred function will depend on the basic sigmoid function, which has the equation $\begin{array}{cc}\sigma \left(x;\alpha ,\beta \right)=\frac{\beta }{1+e^{-\alpha x}},& \text{(4-1)}\end{array}$

[0074] and has a graph of the form shown in FIG. 5.

[0075] The parameter β of the sigmoid function determines a horizontal asymptote for the curve, and the parameter α determines how sharply the function rises as it crosses the y-axis. As α approaches infinity, the sigmoid curve approaches a step function. In the preferred embodiment β=1.0 and α=0.5.

[0076] Because the sigmoid function can transition sharply from a low to a high value, it is a good continuous approximation of discrete processes. Sums of sigmoids can also be used to determine whether or not a variable has a value in an interval. Specifically, for the interval [a, b], define the function

[0077] D(x;a,b)=σ(x−a;α,β)−σ(x−b;α,β)  (4-2).

[0078] Based on the graph of the sigmoid as depicted in FIG. 5, the graph of D(x; a, b) takes the form shown in FIG. 6, which shows the function D(x; a, b) (reference character 60) derived as a sum of two sigmoid functions 62 and 64.

[0079] Since each of the sigmoids 62 and 64 could be made to approximate a step function as closely as desired, the function D(x; a, b) can be defined to very sharply discriminate when x is in the interval [a, b], and can be made to approach a pulse of width b−a as closely as desired.

[0080] Also, since the function D(x; a, b) (reference character 60) approaches zero as x becomes more distant from the interval [a, b], it is possible to sum such interval discriminators (for non-overlapping intervals) and thereby obtain a function which takes a high value when x is in any of the intervals of interest, but is low otherwise. This is shown in FIG. 7 for the two intervals [a₁, b₁] and [a₂, b₂], and it is obvious to those skilled in the art that the construction is generalizable to any finite number of intervals.

[0081] The localizer function 70 shown in FIG. 7 (generated by summing sigmoid functions 72, 74, 76, and 78) can be extended to any finite number of intervals, so such a localizer can be constructed for any corridor of the type in FIG. 1 (one main track, one or more sidings). The sidings are represented along the x-axis between points ai and bi.

[0082] The localizer function 70 has the form $L^{\prime }\left(x;\alpha ,\beta ,a_1,b_1,a_2,b_2\right)=\beta -\sum _{i=1}^2\sigma \left(x-a_i;\alpha ,\beta \right)+\sum _{i=1}^2\sigma \left(x-b_i;\alpha ,\beta \right)$

[0083] The cost function for the scheduling problem of FIG. 2 will be derived below using the localizer function concept, and assuming n_Ssidings. In the preferred embodiment, the cost function will be low if and only if the y-coordinate y_ijfor an intersection of train trajectories lies within the range of a siding, but the localizer function 70 of FIG. 7 in fact displays the opposite effect. Thus we will first definethe localizer function $\begin{array}{cc}{L}^{\prime }\left(x;\alpha ,\beta ,a_1,b_1,\dots ,a_{n_S},b_{n_S}\right)=\beta -\sum _{i=1}^{n_S}\sigma \left(x-a_i;\alpha ,\beta \right)+\sum _{i=1}^{n_S}\sigma \left(x-b_i;\alpha ,\beta \right),& \text{(4-3)}\end{array}$

[0084] which has the desired property of taking a low value if and only if x is in one of the intervals

[0085] [a₁,b₁], . . . ,└a_n,b_n┘,

[0086] and a high value otherwise. That is, Equation 4-3 defines a localizer function that is the inverse of the localizer function 70. See the localizer function 80 in FIG. 8.

[0087] The localizer function as defined above in Equation 4-3 (and taking the form of the inverse of the localizer function 70 in FIG. 7) will now be used to define a cost function which takes lower values as the intersection points of train trajectories are moved toward sidings. Two versions of the cost function are described separately below.

[0088] A Simplified Feasible-Schedule Cost Function

[0089] Now letting n_Tbe the set of all trains to be run in the corridor, and letting L_irepresent the train trajectory line for train T_i(as in FIG. 2). Define a set I of all possible y-coordinates of the intersection points between the train trajectories by

[0090] I={y_ij/{y_ij}=L_i∩L_j& i,jε{1, . . . ,n_T}}  (5-1).

[0091] Note that, with reference to FIG. 2, this set includes all possible intersection points between train trajectories, even though some of those points may not be within the corridor 8 and/or time window of interest. It is necessary to consider such out-of-corridor intersection points because the search process will move the train trajectories, and may bring into the corridor 8 an intersection point that initially was outside the corridor 8.

[0092] To create a cost function that takes on a low value if and only if all intersection points lie within one of the sidings 20, we sum localizer function values derived from Equation 4-3. Specifically define the vector that represents all intersection points in the vicinity of

[0093] {right arrow over (y)}=(y_ij/y_ijεI)  (5-2),

[0094] and define the cost function C′({right arrow over (y)})by $\begin{array}{cc}{C}^{\prime }\left(\stackrel{->}{y}\right)=\sum _{y_{\mathrm{ij}}\in I}{L}^{\prime }\left(y_{\mathrm{ij}};\alpha ,\beta ,a_1,b_1,\dots ,a_{n_S},b_{n_S}\right).& \text{(5-3)}\end{array}$

[0095] The cost function is a multidimensional function of the vector y, where each value of the vector yields a different sum based on localizer function values. Each localizer function value comprising the sum indicates whether an intersection point is in a feasible range (the siding bars 20 of FIG. 2) or not. See the cost function 80 of FIG. 8, where the x-axis represents distance along the corridor. If all of the intersection points relative to a specific siding are feasible, C′({right arrow over (y)}) should take a low value in the vicinity of points represented by that siding; otherwise, it takes a value near the value of β. When many intersection points are involved, β may have to be chosen so that the near-zero sums of a large number of feasible intersection points do not result in a value in the range of β, which would mask the feasibility that is to be discriminated by the function.

[0096] C′({right arrow over (y)}) is a differentiable function of the vector {right arrow over (y)} (the intersection points) and therefore in each of the variables that determine the various intersection points, i.e., departure times and/or speeds of the trains. Therefore the cost function can be used with gradient search technique or other search techniques based on partial derivatives, to minimize the cost function value at sidings. One such technique will be discussed below. Since each intersection point occurring as a component of {right arrow over (y)} is a function of the train departure times and the speeds of the corresponding trains, we may treat the cost function as one which may be optimized by adjusting either speeds or origination times of the trains, or both.

[0097] Accounting for Corridor Endpoints

[0098] The fact that the intersection points in I may not always represent intersections of trajectories within the corridor 8 poses a difficulty for the cost function as defined in Equation 5-3, which is that any intersection point outside the corridor is a “don't care” point for the search process (so long as it remains outside the corridor), but the cost function as defined in Equation 5-3 will assign a high value to such a point. Recall that the cost function of Equation 5-3 is based on the localizer function of Equation 4-3, which is illustrated by reference character 90 in FIG. 9A. Thus as Equation 5-3 is currently formulated, an otherwise feasible solution might be masked by such a “don't care” point.

[0099] In another embodiment of the present invention, the solution involves modifying the localizer function 90. FIG. 9A depicts the localizer function 90, as defined by Equation 4-3, and a modified localizer function 92, which is generated by adding two more sigmoid functions 94 and 96 to account for the end points of the corridor 8. Specifically, define

[0100] e=y-coordinate of the eastern end of the corridor 8,

[0101] w=y-coordinate of the western end of the corridor 8,  (5-4)

[0102] then alter the definition of the localizer function by including the sigmoid functions 74 and 96 as follows.

[0103] L(x;α,β,a₁,b₁, . . . ,a_n,b_n,e,w)=L′(x;α,β,a₁,b₁, . . . ,a_n,b_n)−σ(w−x;α,β)−σ(x−e;α,β)  (5-5).

[0104] The use of the localizer function L of Equation 5-5 also requires rewriting of the cost function in Equation 5-3, as follows. $\begin{array}{cc}C\left(\stackrel{->}{y}\right)=\sum _{y_{\mathrm{ij}}\in I}L\left(y_{\mathrm{ij}};\alpha ,\beta ,a_1,b_1,\dots ,a_{n_S},b_{n_S},e,w\right).& \text{(5-6)}\end{array}$

[0105] This cost function then should take a high value so long as any train trajectory intersection point within the corridor is infeasible, but has a low value for all feasible intersection points, as well as intersection points that fall outside of the corridor.

[0106] Like C′({right arrow over (y)}), C({right arrow over (y)}) is a differentiable function in each component of the vector ({right arrow over (y)}). Any Gradient search techniques or the use of other information based on partial derivatives, can be used to minimize the value of C({right arrow over (y)}) in the regions of the sidings.

[0107] A Balanced Feasible-Schedule Cost Function

[0108] As can be seen from the localizer functions 70, 90, or 92, the sidings (as represented by the x-axis values a₁to b₁) are shown as being of different lengths. In fact, rail corridors typically have sidings of different lengths. The consequence of different length sidings, with respect to the cost function (see Equation (5-6)) is that the cost function minima corresponding to sidings do not have the same y value. See the cost function 80 of FIG. 8. For siding S1, the cost function y value is represented by reference character 82 and the y value for siding S2 is represented by reference character 84. Note that the minimum at reference character 82 has a larger value than the minimum at reference character 84. Because the sidings are different lengths, the sigmoid sum that is creating the minimum uses a narrower portion of the sigmoid function for narrower sidings and, therefore, the associated minimum does not drop down as far as for a wider siding. This effect may cause the cost function gradient optimization process to favor a long siding with a deeper minimum when it is located very close to a short siding with a shallow minimum. In the embodiment discussed below, the cost function will be adjusted to achieve equal minima for all sidings.

[0109] If the derivative of the localizer function 80 has a zero exactly at the midpoint between sidings, then the search process will have no tendency to favor one siding over another. We will call such a localizer function balanced. The situation depicted in FIG. 8 does not assure that the cost function derivative will have a zero properly situated; although the derivative may appear to be zero between the sidings, it can be shown by those skilled in the art, through equation manipulation, that the zero is usually off-center. FIG. 10 illustrates a means to achieve a close approximation to a balanced cost function. In FIG. 10, the intervals [a₁, b₁], [a₂, b₂], and [a₃, b₃] represent locations of sidings along the main corridor. We would like to assure that the derivative of the localizer function, as defined for this corridor, will be zero at the midpoints m₁₂, and m₂₃between sidings. The localizer function generating the cost function is a sum of sigmoids, each of which contributes substantially only within the immediate vicinity of the sidings for which it creates a minimum in the localizer function. If we assume that the localizer function at point m₁₂does not depend significantly on the sigmoid terms other than those used to create minima for the two immediately surrounding sidings, then we may write a simplified localizer in the form

[0110] {tilde over (L)}(x;a₁,b₁,a₂,b₂)=β−σ(x−b₁;α,β)+σ(x−a₁;α,β)−σ(x−b₂;α,β)+σ(x−a₂;α,β)  (5-7).

[0111] Note here that the sigmoid functions used to generate the localizer function are only those sigmoid functions representing sidings to the left and right of the point of interest on the localizer function.

[0112] It can be shown by computation that ${\frac{\partial }{\partial x}\left(\tilde{L}\left(x;a_1,b_1,a_2,b_2\right)\right)|}_{x=m_{12}}=0$

[0113] provided that

[0114] m₁₂−b₁=a₂−m₁₂& m₁₂−a₁=b₂−m₁₂,

[0115] as is shown in FIG. 10. This requirement will force both sidings to be of the same length, of course, and additionally, the next siding corresponding to the interval [a₃, b₃] must then be the same length as the siding corresponding to interval [a₂, b₂]. It follows by induction that all sidings along the corridor must have equal lengths for the localizer function for the corridor to be balanced.

[0116] Bringing such an artifact to bear would have two effects:

[0117] (1) the search might, at least slightly, mislocate intersection points, since the exact position of the sidings would not be reflected in the model;

[0118] (2) siding lengths would not be accurately represented relative to train lengths.

[0119] Of these two drawbacks, the latter is in fact of no consequence, because the modification to the localizer function to account for siding lengths will not affect the subsequent step of the present invention (to be discussed below) wherein train lengths are considered relative to siding lengths. The former effect will be of minor consequence, since getting train intersection points nearly into the vicinity of sidings will allow minor adjustments to train speed to ensure intersections occur at sidings. This step of the present invention will also be discussed further below.

[0120] In another embodiment especially favorable if there is a large discrepancy between the shortest and longest siding, begin with all sidings assumed equal, thereby preventing bias among sidings in the early part of the search, and then adjust the localizer slowly back toward correct siding lengths as the search process iterates. Specifically, this may be implemented in another embodiment of the present invention as follows. Before the search process begins,

[0121] (1) compute the average siding length s_avgas $\begin{array}{cc}{s}_{\mathrm{avg}}=\frac{\sum _{i=1}^{n_S}\left(b_i-a_i\right)}{n_s};& \text{(5-8)}\end{array}$

[0122] (2) redefine the position of each siding S_i(corresponding to corridor interval [a_i, b_i]) as corresponding to the interval [a′_i(0),b′_i(0)], where $\begin{array}{cc}{a}_i^{\prime }=\frac{a_i+b_i-s_{\mathrm{avg}}}{2},\mathrm{and}{b}_i^{\prime }=\frac{a_i+b_i+s_{\mathrm{avg}}}{2};& \text{(5-9)}\end{array}$

[0123] (3) define, for any integer n>0, $\begin{array}{cc}{a}_i^{\prime }\left(n\right)=a_i^{\prime }{e}^{-\lambda n}+a_i\left(1-e^{-\lambda n}\right)\mathrm{and}{b}_i^{\prime }\left(n\right)=b_i^{\prime }{e}^{-\lambda n}+b_i\left(1-e^{-\lambda n}\right),& \text{(5-10)}\end{array}$

[0124] where λ is a positive real number. Note then that $a_i^{\prime }\left(0\right)=a_i^{\prime }\mathrm{and}\underset{n->\infty }{\lim }{a}_i^{\prime }\left(n\right)=a_i,\mathrm{and}$ $b_i^{\prime }\left(0\right)=b_i^{\prime }\mathrm{and}\underset{n->\infty }{\lim }{b}_i^{\prime }\left(n\right)=b_i.$

[0125] Begin the process by letting n=0 and then as the search proceeds, increase n according to some scheme.

[0126] For example, one preferred scheme would be to note when successive values of the cost function (during the gradient search process as discussed below) have a difference smaller than a predetermined threshold (see for example, the threshold value e referred to in conjunction with Equation 8-3 and the textual material following immediately thereafter), then begin to increase n (relative to the differences in siding lengths) and recompute the localizer function until siding lengths are within 5% of being accurate. This will permit the initial localizer to correspond to the balanced localizer, so that sidings will not tend to be favored solely by length. The initial “push” of intersections toward one or another siding will be unbalanced. As n increases, the localizer function will more accurately reflect the true corridor structure, so that eventually an accurate schedule is obtained.

[0127] Accounting for Train Lengths Versus Siding Lengths

[0128] The cost function as described above permits a search for a feasible schedule only insofar that trains will meet in the vicinity of sidings. No reference has been made to the lengths of the trains relative to the sidings, and if two trains have a “feasible” meet at a siding that will hold neither of them, then the situation is not actually feasible. There are other reasons that trains may not use a siding, related to grade, transportation of hazardous materials, etc., so the following analysis to block the use of a siding by a given train refers to more situations than just train length versus siding length.

[0129] The cost function of Equation 5-6 will not prevent such an infeasibility from occurring, but in another embodiment, a simple modification of the localizer functions (Equation 5-5) on which the cost function is based will suffice to prevent such infeasibilities.

[0130] In particular, the cost function contains a term for each possible train trajectory intersection point. In the previous embodiment all such terms are of exactly the same form. Now suppose that we define the localizer functions to be specific to each possible intersection point of train trajectories, as follows. In this case, we generalize from the context of FIG. 2, and assume a total of n_Ssidings S₁, . . . ,S_nalong the corridor, and n_Ttrains T₁, . . . , T_n. We need the following notation:

[0131] let

[0132] H_i=the length of siding S_i(i=1, . . . , n_S)  (6-1),

[0133] and let

[0134] M_i=the length of train T_i(i=1, . . . ,n_T)  (6-2).

[0135] For any two trains T_iand T_j, define the following set of sidings from amongst all sidings in the corridor:

[0136] S_ij={S_k/kε{1, . . . , N} & ((M_i≦H_k)v(M_j≦H_k))}  (6-3).

[0137] S_ijis the subset of sidings along the corridor on which at least one of the two trains T_iand T_jcan be sided. Now, if the localizer function for the intersection point of the train trajectories of T_iand T_jdoes not include the sigmoid terms (see Equation 5-7) corresponding to sidings not in S_ij, then it will remain high even though y_ijis within a siding, but the siding is too short for either train. In this way, accounting for siding versus train length actually reduces the computational complexity of the cost function.

[0138] To specifically redefine the cost function in this form, first redefine the localizer functions to be specific to train pairs, i.e., $\begin{array}{cc}{L}_{\mathrm{ij}}\left(y_{\mathrm{ij}};\alpha ,\beta \right)=\beta -\sum _{h\in S_{\mathrm{ij}}}\sigma \left(y_{\mathrm{ij}}-a_h;\alpha ,\beta \right)+\sum _{h\in S_{\mathrm{ij}}}\sigma \left(y_{\mathrm{ij}}-b_h;\alpha ,\beta \right)-\sigma \left(w-y_{\mathrm{ij}};\alpha ,\beta \right)-\sigma \left(y_{\mathrm{ij}}-e;\alpha ,\beta \right).& \text{(6-4)}\end{array}$

[0139] Where the subscript “h” identifies a siding.

[0140] Finally redefine the cost function as $\begin{array}{cc}C\left(\stackrel{->}{y}\right)=\sum _{y_{\mathrm{ij}}\in I}{L}_{\mathrm{ij}}\left(y_{\mathrm{ij}};\alpha ,\beta \right),& \text{(6-6)}\end{array}$

[0141] which extends the definition of feasibility so that now the value of C({right arrow over (y)}) will be low if and only if

[0142] (1) all train trajectory intersections occur on siding bars, and

[0143] (2) at least one of the two trains in such an intersection can be sided in the corresponding siding.

[0144] Note that this technique can be extended beyond the consideration of train length versus siding length: if neither of two trains T_iand T_jcan be sided on siding S_kfor any reason, then the localizer for the intersection point y_ijshould omit the term corresponding to S_k. For example, we may have a case where a coal train could be sided at S_k, but would be unable to restart because of grade, but the interfering train, a multimodal, is absolutely not to be sided for a coal train. In this case, the siding may be long enough for either train, but would be precluded from consideration anyway. Clearly in other embodiments the definition of each S_ijcan be contracted to exclude cases such as this, thereby sharpening the ability of the search process to prevent unacceptable sidings.

[0145] Economic Costs, Early Departure, and Speed Constraints

[0146] The cost function as described by either Equation 5-6 or 6-6 will facilitate the finding of feasible train schedules, but includes no cognizance of the other effects of altering individual train schedules to achieve feasibility. In another embodiment, the cost function is modified so that it jointly considers schedule feasibility, and the economic cost of late arrival.

[0147] Economic Costs (i.e., Late Arrival) Function

[0148] Railroad freight service may incur various types of incentives for on-time delivery of freight. For the moment, consider just two types of delay penalties:

[0149] (1) step function penalty—if a train T_imisses a preset delivery time t_i, there is a fixed penalty cost h_i;

[0150] (2) step function plus linear increase—if the preset delivery time t, is missed, there is an immediate penalty h_i(possibly 0) which thereafter linearly increases at a rate of m_idollars per hour.

[0151] FIG. 11A depicts a single generic form for both of these cases, since both h_iand m_imay be zero or positive. Thus, FIG. 11A illustrates a combined penalty function including both a step penalty plus a linear penalty.

[0152] The cost function as proposed is in not a differentiable function since it lacks a defined slope at the time t_i. This fact precludes, or at least complicates, use of any gradient search technique for minimizing economic cost unless special allowances are made at or near the time t_i. For this reason, FIGS. 11B and 11C depict two approximations to the cost function, a step plus linear penalty, and a linear penalty only, respectively. In both figures, a line segment is grafted onto a sigmoid function in such a way that the resulting function remains differentiable at all points.

[0153] For the step plus linear penalty, a sigmoid is used to represent the cost up to a time slightly beyond t_i, to which is then appended a line of slope m_i. See FIG. 11B. Provided the crossover point from sigmoid to line segment is chosen at the point of the sigmoid where the slope is exactly m_i, (reference character 110) the resulting approximation is differentiable at all points, and therefore smoothly integratable into a gradient search process. If (t_c, y_c) represents the crossover point; then the differentiable version of the penalty function may be defined by $\begin{array}{cc}{A}_i\left(t;t_i,h_i,m_i\right)=\{\begin{array}{cc}{h}_i\sigma \left(t-t_i;{\alpha }_i\right)& \mathrm{for}t<t_c\\ m_{i}t+y_c-m_i{t}_i& \mathrm{for}t\ge t_c\end{array}& \text{(7-1)}\end{array}$

[0154] The sigmoids used here will all have β_lvalues of 1, so that the notation for the parameter β_lin each sigmoid will be suppressed. This choice is made so that the asymptote of the sigmoid is determined to be h_i, in conjunction with the penalty value to be represented.

[0155] The value of α_iis positive, and may be chosen to approximate the step cost as sharply as desired. In one embodiment the search is started with “gentle” sigmoids, then increase the values of the α_i's as the search progresses. This allows the early search to progress toward correct economic decisions rapidly, and then in the later stages of search, the information concerning economic cost is sharpened to provide more accurate final results.

[0156] In order to determine the crossover point 110 ((t_c, y_c) in FIG. 11B), it is necessary to solve the equation $\begin{array}{cc}{\frac{\partial }{\partial t}\left(h_i\sigma \left(t;{\alpha }_i\right)\right)|}_{t=t_c}=m_i& \text{(7-2)}\end{array}$

[0157] for the value of t_c, with t_c>t_i. The technique for solving this equation is well known to those skilled in the art. It should be mentioned that the slope of σ(t−t_i;α_i) is everywhere positive, and takes a maximum at the point t=t_i. That maximum can be driven as high as possible by selecting a large α_i, so solving Equation 7-2 is always possible.

[0158] Finally, for the purposes of expressing the gradient as will be explained below, note that the independent variable t in Equation 7-1 is in fact a function of the departure time d_iand speed s_iof train T_i, and therefore we may rewrite the equation as $\begin{array}{cc}{A}_i\left(s_i,d_i;t_i,h_i,m_i\right)=\{\begin{array}{cc}{h}_i\sigma \left(d_i+\frac{L}{s_i}-t_i;{\alpha }_i\right)& \mathrm{for}{d}_i+\frac{L}{s_i}<t_c\\ m_i{d}_i+\frac{L}{s_i}+y_c-m_i{t}_i& \mathrm{for}{d}_i+\frac{L}{s_i}\ge t_c.\end{array}& \text{(7-3)}\end{array}$

[0159] FIG. 10C also uses a transition from sigmoid to line segment at point 112 on the sigmoid where the slope is exactly that of the line: the difference is that in this case the crossover point t_cis less than t_i. Except for that fact, the approximating function has a description identical to that provided in Equations 7-1 and 7-3.

[0160] Now, we extend the cost function of Equation 5-6 or Equation 6-6 as follows. The extended cost function accounting for both schedule feasibility and economic cost is defined by $\begin{array}{cc}F\left(\overrightarrow{y}\right)=\eta C\left(\overrightarrow{y}\right)+\left(1-\eta \right)\sum _{i=1}^{n_T}{A}_i\left(\overrightarrow{s},\overrightarrow{d};t_i,h_i,m_i\right),& \text{(7-4)}\end{array}$

[0161] where

[0162] ηε[0,1] is a weighting factor between 0 and 1 used to adjust the relative importance between economic and schedule feasibility considerations.

[0163] {right arrow over (d)}=(d₁,d₂, . . . ,d_n) is the vector of train departure times,

[0164] {right arrow over (s)}=(s₁,s₂, . . . ,s_n) is the vector of train speeds.

[0165] In fact, the intersection points y of train trajectories are functions of the train departure times and speeds, so we may rewrite Equation 7-4 in the form $\begin{array}{cc}F\left(\overrightarrow{s},\overrightarrow{d}\right)=\eta C\left(\overrightarrow{s},\overrightarrow{d}\right)+\left(1-\eta \right)\sum _{i=1}^{n_T}{A}_i\left(s_i,d_i;t_i,h_i,m_i\right),& \text{(7-5)}\end{array}$

[0166] and it is from this latter form that the gradient may be directly computed as discussed below.

[0167] The value of the weighting factor η must be chosen, and the choice is of some importance. Note that the cost function as defined in Equation 7-5 will be driven upward both by infeasible scheduling choices, as well as choices that make trains late, and vice versa. The difficulty arises when changes of departure times or speeds cause countervailing effects in the two halves of the cost function of Equation 7-5. If the first term, representing feasibility, is driven up by less than the second term, representing timeliness, is driven down, then the search process may be emphasizing economic cost to such an extent that it converges on infeasible schedules.

[0168] In one embodiment, the weighting actor η can be varied during the search. For example, starting with a low value of η would tend to try to force low economic cost at the expense of feasibility. This might cause the trains to swap places in the lineup, to improve the overall timeliness of arrivals, before the actual emphasis begins on selecting speeds and departure times that create a feasible schedule. In any event, the decision as to how to vary η during the search will benefit from actual testing with examples, and fmal mechanism for modulating η will necessarily come from experience familiar to those skilled in the art.

[0169] An approximate process for gauging the weighting factor η is to note that the cost components C({right arrow over (s)},{right arrow over (d)}) and $\sum _{i=1}^{n_T}A\left(s_i,d_i;t_i,h_i,m_i\right)$

[0170] comprise different numbers of summands, and therefore have different magnitudes approximately in proportion to the number of summands involved. For example, if there are a total of twenty trains, resulting in sixty intersections on the string graph, then C({right arrow over (s)},{right arrow over (d)}) comprises sixty summands and $\sum _{i=1}^{n_T}A\left(s_i,d_i;t_i,h_i,m_i\right)$

[0171] comprises twenty summands. To more or less equalize the effects of these two contributions to the cost function, one would set the weight η to the value η=20/(60+20)=0.25, thereby equalizing the contribution of each half of the cost function (i.e., the two terms C({right arrow over (s)},{right arrow over (d)}) and $\sum _{i=1}^{n_T}A\left(s_i,d_i;t_i,h_i,m_i\right))$

[0172] to the total cost. From this example, it can be seen that the establishment of a specific value of η is very specific to the situation under study, as is generally recognized by those proficient in the art of complex optimization.

[0173] Early Departure Cost Function

[0174] The late penalty assessed for economic reasons will tend to prevent train departures from being arbitrarily late. However, the formulations of cost functions so far given (Equations 5-6, 6-6, 7-5) have no terms which prevent the train trajectories from being arbitrarily early. A cost function to prevent early departures can be formulated in terms of the ubiquitous sigmoid function by defining a cost $\begin{array}{cc}E\left(\overrightarrow{d}\right)=\sum _{i=1}^{n_T}\left[1-\sigma \left(d_i-e_i;{\alpha }_i^{\prime }\right)\right],& \text{(7-6)}\end{array}$

[0175] where

[0176] e_iis the earliest possible departure time for train T_i, and

[0177] α′_iis denoted with a prime to distinguish it from the α_iof Equation 7-3.

[0178] FIG. 12 represents a term of this cost function for train T_i; clearly, it rapidly becomes high as train T_iis pushed toward an unrealizable departure time, and rapidly drops as the departure time enters the realizable region. There is no actual economic cost associated with early departures, just a feasibility issue. Therefore the terms of Equation 7-6 representing each train are arbitrarily given a height of 1, (i.e., sigmoid asymptotic value of 1) and this Equation 7-6 can likewise be combined with the cost functions for schedule feasibility and economic cost. Specifically, let

[0179] G({right arrow over (s)},{right arrow over (d)})=η₁C({right arrow over (s)},{right arrow over (d)})+η₂A({right arrow over (s)},{right arrow over (d)})+η₃E({right arrow over (d)})  (7-7),

[0180] where

[0181] η₁+η₂+η₃=1  (7-8).

[0182] In one embodiment, the specific weighting of the components of cost in Equation 7-7 can be calculated as discussed above in the example with twenty trains and sixty intersections points in the corridor. The schedule feasibility and early departure terms will each have sixty summands and the economic penalty term will have twenty terms. Using an equation similar to the one set forth above for calculating η, we calculate η₁=1/7, η₂=3/7 and η₃=3/7. Other weighting values can be established based on specific user circumstances.

[0183] Maximum Train Speed Cost Function

[0184] In the embodiment when the search process is permitted to vary train speeds in order to achieve feasibility and cost minimization, there must be a means to prevent the speeds from exceeding practical limits for the trains and tracks involved. In this embodiment we will create an additional component of the cost function that will enforce such speed constraints. Such a speed constraint can be implemented analogously to the early departure constraint of Equation 7-6. Specifically, define a speed cost function as $\begin{array}{cc}V\left(\overrightarrow{s}\right)=\sum _{i=1}^{n_T}\sigma \left(s_i-s_i^{\left(\max \right)}\right),& \text{(7-9)}\end{array}$

[0185] where

[0186] s_i^(max)=the maximum allowable speed for train T_i.

[0187] Like the other cost functions discussed herein, since the maximum speed cost function is derived from a sum of sigmoid functions, it id a differentiable function with respect to the intersection points of the trains on the corridor. Therefore a gradient search process can be used to find the minima of the cost function values. The total cost function, including feasibility of meets and passes, constraints on early departures and late arrivals (i.e., economic penalty), as well as constraints on maximum train speed then is a generalization of Equation 7-7, namely

[0188] G({right arrow over (s)},{right arrow over (d)})=η₁C({right arrow over (s)},{right arrow over (d)})+η₂A({right arrow over (s)},{right arrow over (d)})+η₃E({right arrow over (d)})+η₄V({right arrow over (s)})  (7-10),

[0189] where

[0190] η₁+η₂+η₃+η₄=1  (7-11).

[0191] The specific values of the weighting factors for the components of Equation 7-10 can be determined by experiment. In one embodiment, using the same scheme as set forth above in conjunction with Equation (7-8), for twenty trains and sixty intersections, η₁=0.1 and η₂=η₃=η₄=0.3.

[0192] The Gradient Search Process

[0193] The gradient Δ𝑓({right arrow over (x)}) of any function 𝑓({right arrow over (x)}) is a vector in the same space as the independent variable x which points in direction of maximum change of 𝑓({right arrow over (x)}) within a small local area on the function's surface, thereby pointing the way toward a local minimum or maximum. As such, it is much heralded in the legends and poetry of optimization theory. Calculation of the gradient of the various cost functions discussed below will permit location of the local minima identifying schedule feasibility. In the current context of train scheduling, as will be appreciated by those skilled in the art, there are a number of possible parameters describing a train trajectory that may be varied to resolve conflicts within a rail corridor, i.e., to drive the cost function lower. The mathematics for a gradient search varying only the departure times or speeds of trains, and then varying both departure times and train speeds is discussed below. We will first deal only with the cost function associated with schedule feasibility (Equation 5-6) but will then extend the cost function to include considerations of economic costs, early departures and maximum train speed as discussed above, and represented by the cost function of Equation 7-10.

[0194] Gradient Search to Optimize Schedule Feasibility by Varying Only the Train Departure Times

[0195] First, assume that there are n_Ttrains and permit the vector {right arrow over (y)} (as represented in Equation (5-2)) to contain all possible intersection points. But each intersection point y_ijhas the characterization given in Equation 3-6, which is repeated here for convenience. $\begin{array}{cc}{y}_{\mathrm{ij}}=\frac{s_i{s}_j\left(d_i-d_j\right)+\left(s_i{\theta }_j-s_j{\theta }_i\right)L}{s_i-s_j},& \text{(8-1)}\end{array}$

[0196] then y_ijis directly expressed in terms of departure times and speeds for all of the trains in the schedule. Also for convenience, recall the notation for speeds and departure times originally introduced above, which are repeated below.

[0197] L=the length of the corridor,

[0198] s_i=the speed of train T_i(taken as a negative value for T_iwestbound),

[0199] d_i=the departure time (time of entry into corridor) of train t_i, and ${\theta }_i=\{\begin{array}{cc}0& \mathrm{for}{T}_i\mathrm{eastbound}\\ 1& \mathrm{for}{T}_i\mathrm{westbound}.\end{array}$

[0200] Next define the vectors {right arrow over (s)}=(s₁, . . . ,s_n) and {right arrow over (d)}=(d₁, . . . ,d_n). Express the cost function in the following terms. For notational convenience, suppress the dependence of the localizer and cost functions on α and β. $\begin{array}{cc}C\left(\overrightarrow{y}\right)=C\left(\overrightarrow{s},\overrightarrow{d}\right)=\sum _{y_{\mathrm{ij}}\in I}{L}_{\mathrm{ij}}\left(y_{\mathrm{ij}}\right)=\sum _{y_{\mathrm{ij}}\in I}{L}_{\mathrm{ij}}\left(s_i,s_j,d_i,d_j\right).& \text{(8-2)}\end{array}$

[0201] The objective is to vary the vector {right arrow over (d)} (train departure times) in order to drive the cost function lower, and one technique which can at least locate a local minimum of the cost function is the gradient-directed descent, defined iteratively as follows.

[0202] (1) Start with an initial estimate for departure time, {right arrow over (d)}₀for each train n_T, a stopping criterion, ε>0, and a step size h.

[0203] (2) For the estimate, {right arrow over (d)}_n, compute the gradient ∇^{({right arrow over (d)})}C({right arrow over (d)})|_{{right arrow over (d)}={right arrow over (d)}}of the cost function at {right arrow over (d)}_n, for varying only {right arrow over (d)}, and normalize it so that it has an absolute value of 1, i.e., define $\begin{array}{cc}\overrightarrow{g}=\frac{{{\nabla }^{\left(\overrightarrow{d}\right)}C\left(\overrightarrow{d}\right)|}_{\overrightarrow{d}={\overrightarrow{d}}_n}}{\left|{\nabla }^{\left(\overrightarrow{d}\right)}C\left(\overrightarrow{d}\right){❘}_{\overrightarrow{d}={\overrightarrow{d}}_n}\right|}.& \text{(8-3)}\end{array}$

[0204] In the notation, the dependence of the cost function on {right arrow over (s)} is suppressed, since for the nonce we are varying only {right arrow over (d)}.

[0205] (3) Compute the value C_n=C({right arrow over (d)}_n), compute {right arrow over (d)}_n+1={right arrow over (d)}_n−h{right arrow over (g)}, and then compute C_n+1=C({right arrow over (d)}_n+1).

[0206] (4) If |C_n−C_n+1|<ε, then the search is stopped, and {right arrow over (d)}_n+1is accepted as the final answer. Otherwise, replace {right arrow over (d)}_nwith {right arrow over (d)}_n+1and return to Step (2). In the preferred embodiment, the search is stopped when |C_n−C_n+1|≦(0.001)|C₀−C₁. The stopping threshold for such problems is very situation dependent, as is generally recognized by practitioners of the art of optimization.

[0207] It remains to explicitly represent the gradient ∇^{({right arrow over (d)})}C({right arrow over (d)})|_{{right arrow over (d)}={right arrow over (d)}}which is used in the iteration. The cost function as shown in Equation 5-6 is a function of the vector of intersection points, {right arrow over (y)}, and the components of {right arrow over (y)} are functions of the components of the vectors {right arrow over (s)} and {right arrow over (d)}. Since at this point only {right arrow over (d)} is variable, the gradient ∇^{({right arrow over (d)})}C({right arrow over (d)}) of the cost function, is a vector of the form $\begin{array}{cc}{\nabla }^{\left(\overrightarrow{d}\right)}C\left(\overrightarrow{d}\right)=\left(\frac{\partial }{\partial d_1}C\left(\overrightarrow{d}\right),\frac{\partial }{\partial d_2}C\left(\overrightarrow{d}\right),\dots ,\frac{\partial }{\partial d_{n_T}}C\left(\overrightarrow{d}\right)\right),& \text{(8-4)}\end{array}$

[0208] and we may obtain each component of ∇^{({right arrow over (d)})}C({right arrow over (d)}) by applying the chain rule for differentiation: $\begin{array}{cc}\begin{array}{cc}\frac{\partial }{\partial d_k}C\left(\overrightarrow{d}\right)=\frac{\partial }{\partial d_k}\left(\sum _{y_{\mathrm{ij}}\in I}{L}_{\mathrm{ij}}\left(y_{\mathrm{ij}}\right)\right)=\frac{\partial }{\partial d_k}\left(\sum _{y_{\mathrm{ij}}\in Ii=k⩔j=k}{L}_{\mathrm{ij}}\left(y_{\mathrm{ij}}\right)\right)& \begin{array}{c}\text{(8-5)}\\ \end{array}\\ =\sum _{y_{\mathrm{ik}}\in I}\left[\frac{\partial }{\partial y_{\mathrm{ik}}}\left(L_{\mathrm{ik}}\left(y_{\mathrm{ik}}\right)\right)\frac{\partial }{\partial d_k}\left(y_{\mathrm{ik}}\right)\right]& \\ =\sum _{y_{\mathrm{ik}}\in I}\frac{\partial }{\partial y_{\mathrm{ik}}}\left(L_{\mathrm{ik}}\left(y_{\mathrm{ik}}\right)\right)\frac{\partial }{\partial d_k}\left(\frac{\begin{array}{c}{s}_i{s}_k\left(d_i-d_k\right)+\\ \left(s_i{\theta }_k-s_k{\theta }_i\right)L\end{array}}{s_i-s_k}\right)& \\ =\sum _{y_{\mathrm{ik}}\in I}\frac{\partial }{\partial y_{\mathrm{ik}}}\left(L_{\mathrm{ik}}\left(y_{\mathrm{ik}}\right)\right)\left(\frac{-s_i{s}_k}{s_i-s_k}\right).& \begin{array}{c}\text{(8-6)}\\ \end{array}\end{array}& \end{array}$

[0209] Using Equation 8-6 and Lemma A3 in the Appendix A we can finally express the k-th component of the gradient as $\begin{array}{cc}\begin{array}{c}\frac{\partial }{\partial d_k}\left(C\left(\overrightarrow{d}\right)\right)=\frac{\alpha }{\beta }\sum _{y_{\mathrm{ik}}\in I}\{\sum _{j=1}^{n_S}[\sigma \left(y_{\mathrm{ik}}-b_j\right)\left(\beta -\sigma \left(y_{\mathrm{ik}}{b}_j\right)\right)-\\ \sigma \left(y_{\mathrm{ik}}-a_j\right)\left(\beta -\sigma \left(y_{\mathrm{ik}}-a_j\right)\right)+\\ \sigma \left(w-y_{\mathrm{ik}}\right)\left(\beta -\sigma \left(w-y_{\mathrm{ik}}\right)\right)-\\ \sigma \left(y_{\mathrm{ik}}-e\right)\left(\beta -\sigma \left(y_{\mathrm{ik}}-e\right)\right)]\left(\frac{s_i{s}_k}{s_i-s_k}\right)\}.\end{array}& \text{(8-7)}\end{array}$

[0210] Gradient Search to Optimize Schedule Feasibility by Varying Only the Train Speeds

[0211] Much of what was developed above can be applied here as well. The primary difference is that we now emphasize that C({right arrow over (y)}) may be regarded as a function of the vector {right arrow over (s)}, with {right arrow over (d)} held constant, and we wish to vary {right arrow over (s)} to seek a local minimum of the cost function, and suppress the dependence on {right arrow over (d)}. We may therefore represent C({right arrow over (y)}) as

[0212] C({right arrow over (y)})=C({right arrow over (s)})  (8-8).

[0213] Starting with $\begin{array}{cc}{\nabla }^{\left(\stackrel{->}{s}\right)}C\left(\stackrel{->}{s}\right)=\left(\frac{\partial }{\partial s_1}C\left(\stackrel{->}{s}\right),\frac{\partial }{\partial s_2}C\left(\stackrel{->}{s}\right),\dots ,\frac{\partial }{\partial s_{n_T}}C\left(\stackrel{->}{s}\right)\right),& \text{(8-9)}\end{array}$

[0214] we need to compute the components of the form $\frac{\partial }{\partial s_k}\left(C\left(\stackrel{->}{s}\right)\right)$

[0215] using the differentiation chain rule. To that end, note that $\begin{array}{cc}\begin{array}{c}\frac{\partial }{\partial s_k}\left(C\left(\stackrel{->}{s}\right)\right)=\frac{\partial }{\partial s_k}\left(\sum _{y_{\mathrm{ij}}\in I}{L}_{\mathrm{ij}}\left(y_{\mathrm{ij}}\right)\right)=\frac{\partial }{\partial s_k}\left(\sum _{y_{\mathrm{ik}}\in I}{L}_{\mathrm{ik}}\left(y_{\mathrm{ik}}\right)\right)\\ =\sum _{y_{\mathrm{ij}}\in I}\frac{\partial }{\partial y_{\mathrm{ik}}}\left(L_{\mathrm{ik}}\left(y_{\mathrm{ik}}\right)\right)\frac{\partial }{\partial s_k}\left(y_{\mathrm{ik}}\right),\end{array}& \text{(8-10)}\end{array}$

[0216] and we proceed to obtain an explicit expression for $\frac{\partial }{\partial s_k}\left(y_{\mathrm{ik}}\right),$

[0217] as follows. $\begin{array}{cc}\begin{array}{c}\frac{\partial }{\partial s_k}\left(y_{\mathrm{ik}}\right)=\frac{\partial }{\partial s_k}\left(\frac{s_i{s}_k\left(d_i-d_k\right)+\left(s_i{\theta }_k-s_k{\theta }_i\right)L}{s_i-s_k}\right)\\ =\frac{\begin{array}{c}\left[s_i{s}_k\left(d_i-d_k\right)+\left(s_i{\theta }_k-s_k{\theta }_i\right)\right]\left(-1\right)+\\ \left(s_i-s_k\right)\left[s_i\left(d_i-d_k\right)-{\theta }_{i}L\right]\end{array}}{{\left(s_i-s_k\right)}^2}\\ =\frac{s_i\left(d_i-d_k\right)\left(s_i-2s_k\right)+L\left(2s_k{\theta }_i-s_i{\theta }_k-s_i{\theta }_i\right)}{{\left(s_i-s_k\right)}^2}.\end{array}& \text{(8-11)}\end{array}$

[0218] Exploiting Equations 8-10, 8-11, and Lemma A3 of the Appendix A provides a final explicit form for the gradient, as shown below $\begin{array}{cc}\begin{array}{c}\frac{\partial }{\partial s_k}\left(C\left(\stackrel{->}{s}\right)\right)=\frac{\alpha }{\beta }\sum _{y_{\mathrm{ik}}\in I}\{\sum _{j=1}^N[\sigma \left(y_{\mathrm{ik}}-b_j\right)\left(\beta -\sigma \left(y_{\mathrm{ik}}-b_j\right)\right)-\\ \sigma \left(y_{\mathrm{ik}}-a_j\right)\left(\beta -\sigma \left(y_{\mathrm{ik}}-a_j\right)\right)+\\ \sigma \left(w-y_{\mathrm{ik}}\right)\left(\beta -\sigma \left(w-y_{\mathrm{ik}}\right)\right)-\\ \sigma \left(y_{\mathrm{ik}}-e\right)\left(\beta -\sigma \left(y_{\mathrm{ik}}-e\right)\right)]\frac{\partial }{\partial s_k}\left(y_{\mathrm{ik}}\right)\}.\end{array}& \text{(8-12)}\end{array}$

[0219] The search rule using the gradient as computed in Equation 8-12 is an exact analogue of the search rule given in Equation 8-7 with any occurrence of the vectors {right arrow over (d)},{right arrow over (d)}₀,{right arrow over (d)}_n,{right arrow over (d)}_n+1. . . replaced with the vectors {right arrow over (s)},{right arrow over (s)}₀,{right arrow over (s)}_n{right arrow over (s)}_n+1. . . , respectively.

[0220] Gradient Search to Optimize Schedule Feasibility by Varying Both Departure Times and Train Speeds

[0221] Bearing in mind that the speeds and departure times of trains can be varied independently, we may also exploit the expression of the cost function as a function of both {right arrow over (s)} and {right arrow over (d)}, i.e.,

[0222] C({right arrow over (y)})=C({right arrow over (s)},{right arrow over (d)}),

[0223] and consider the joint variation of speed and departure time to seek a cost function local minimum. In this case, the gradient vector takes the form $\begin{array}{cc}{\nabla }^{\left(\stackrel{->}{s},\stackrel{->}{d}\right)}C\left(\stackrel{->}{s},\stackrel{->}{d}\right)=\left(\frac{\partial }{\partial s_1}C\left(\stackrel{->}{s},\stackrel{->}{d}\right),\dots ,\frac{\partial }{\partial s_{n_T}}C\left(\stackrel{->}{s},\stackrel{->}{d}\right),\frac{\partial }{\partial d_1}C\left(\stackrel{->}{s},\stackrel{->}{d}\right),\dots ,\frac{\partial }{\partial d_{n_T}}C\left(\stackrel{->}{s},\stackrel{->}{d}\right)\right).& \text{(8-13)}\end{array}$

[0224] Because {right arrow over (s)} and {right arrow over (d)} are not functionally dependent on each other, it follows that $\begin{array}{cc}\frac{\partial }{\partial s_k}\left(C\left(\stackrel{->}{s},\stackrel{->}{d}\right)\right)=\frac{\partial }{\partial s_k}\left(C\left(\stackrel{->}{s}\right)\right)\mathrm{and}\frac{\partial }{\partial d_k}\left(C\left(\stackrel{->}{s},\stackrel{->}{d}\right)\right)=\frac{\partial }{\partial d_k}\left(C\left(\stackrel{->}{d}\right)\right),& \text{(8-14)}\end{array}$

[0225] so the components of the gradient of Equation 8-13 are aSeady determined by Equations 8-7, 8-1 1, and 8-12.

[0226] The search rule in this case is of the same form as Equation 8-7, except that we consider the aggregate vector

[0227] {right arrow over (v)}=({right arrow over (s)},{right arrow over (d)})  (8-15),

[0228] and replace all references to {right arrow over (d)},{right arrow over (d)}₀,{right arrow over (d)}_n,{right arrow over (d)}_n+1in that rule with references to {right arrow over (v)},{right arrow over (v)}₀,{right arrow over (v)}_n{right arrow over (v)}_n+1, respectively.

[0229] Including Early Departure Effects in the Gradient Search

[0230] Recall from the discussion above that a cost function causing high cost for early train departures, and low cost otherwise, can be posed in terms of the sigmoid function. Repeating Equation 7-6, $\begin{array}{cc}E\left(\stackrel{->}{d}\right)=\sum _{i=1}^{n_T}\left[1-\sigma \left(d_i-e_i;{\alpha }_i^{\prime }\right)\right],& \text{(8-16)}\end{array}$

[0231] where

[0232] e_iis the earliest possible departure time for train T_i, and

[0233] α′_iaffects the steepness of the rise of cost as early departure is approached.

[0234] The value of α′ may be set by experiment, but the results should not be particularly sensitive to its value. A good first guess, in one embodiment, for the value of α′_iwould be 0.8, although this parameter might be made smaller if there is some latitude as to the earliest departure times.

[0235] If we wish to combine this early departure cost with schedule feasibility cost, we do so in a weighted sum of terms, i.e.,

[0236] D({right arrow over (s)},{right arrow over (d)})=ηC({right arrow over (s)},{right arrow over (d)})+(1η)E({right arrow over (d)})ηε(0,1)  (8-17),

[0237] This concept was previously discussed above. See for example, Equation 7-7 where the aggregate cost function includes schedule feasibility, economic cost, and early departure effects.

[0238] Since the gradient operation is linear on the space of functions to which it applies, we may write

[0239] ∇(D({right arrow over (s)},{right arrow over (d)}))=η∇(C({right arrow over (s)},{right arrow over (d)}))+(1−η)∇(E({right arrow over (d)}))  (8-18).

[0240] We will rely on the previous gradient computations for the first term on the right side of Equation 8-17. See Equation 8-7 with the substitutions set forth in Equation 8-15 and the text following.

[0241] To deal with the second term on the right side of Equation 8-17 or 8-18, assume only the departure time vector {right arrow over (d)} will be varied in a search for a schedule which is both feasible and prevents early departures. We then wish to determine the gradient of E({right arrow over (d)}) relative to the vector {right arrow over (d)}, which is of the form $\begin{array}{cc}{\nabla }^{\left(\stackrel{->}{d}\right)}\left(E\left(\stackrel{->}{d}\right)\right)=\left(\frac{\partial }{\partial d_1}\left(E\left(\stackrel{->}{d}\right)\right),\frac{\partial }{\partial d_2}\left(E\left(\stackrel{->}{d}\right)\right),\dots ,\frac{\partial }{\partial d_{n_T}}\left(E\left(\stackrel{->}{d}\right)\right)\right),& \text{(8-19)}\end{array}$

[0242] and (see Equations 8-6 and A-2) $\begin{array}{cc}\begin{array}{c}\frac{\partial }{\partial d_k}\left(E\left(\stackrel{->}{d}\right)\right)=\left(1-\eta \right)\left[\sum _{i=1}^{n_T}\frac{\partial }{\partial d_k}\left(1-\sigma \left(d_i-e_i;{\alpha }_i^{\prime }\right)\right)\right]\\ =-\left(1-\eta \right)\left(1-\sigma \left(d_k-e_k;{\alpha }_k^{\prime }\right)\sigma \left(d_k-e_k;{\alpha }_k^{\prime }\right)\right).\end{array}& \text{(8-20)}\end{array}$

[0243] We can now construct the gradient ∇^{({right arrow over (d)})}(D({right arrow over (s)},{right arrow over (d)})) using Equations 8-7, 8-18, and 8-20. Departure times are independent of train speeds, so the cost component E({right arrow over (d)}) does not depend on the speeds {right arrow over (s)}. Thus the final form of the gradient

[0244] ∇^{({right arrow over (s)},{right arrow over (d)})}E({right arrow over (s)},{right arrow over (d)})=(E₁, . . . ,E_n,E_n+1, . . . ,E_2n)  (8-21)

[0245] with both train speeds and departure times variable, can be summarized as $\begin{array}{cc}{E}_i=\{\begin{array}{c}\eta \frac{\partial }{\partial s_i}\left(C\left(\stackrel{->}{s},\stackrel{->}{d}\right)\right)\mathrm{for}i\in \left\{1,\dots ,n_T\right\}\\ \eta \frac{\partial }{\partial d_i}\left(C\left(\stackrel{->}{s},\stackrel{->}{d}\right)\right)-\left(1-\eta \right)\left(1-\sigma \left(d_i-e_i;{\alpha }_i^{\prime }\right)\sigma \left(d_i-e_i;{\alpha }_i^{\prime }\right)\right)\\ \mathrm{for}i\in \left\{n_T+1,\dots ,2n_T\right\}\end{array}& \text{(8-22)}\end{array}$

[0246] where reference is implicitly made to Equations 8-7, 8-12, and 8-15.

[0247] Including the Economic Costs in the Gradient Search

[0248] The types of costs incurred by railroads for late deliveries were discussed above, and there was provided a differentiable approximation to the function of late costs expressed as a function of time. By using such an approximation, which is everywhere differentiable, the avoidance of late costs can be incorporated into the gradient search process. Arrival times are affected by both train speeds and train departure times, although either speed, departure time, or both may be varied during the search.

[0249] The form of the late cost approximation function was given by (see Equation 7-3) $\begin{array}{cc}{A}_i\left(u_i;t_i,h_i,m_i\right)=\{\begin{array}{c}{h}_i\sigma \left(u_i-t_i;{\alpha }_i\right)\mathrm{for}{u}_i<t_c\\ y_c+m_i\left(u_i-t_i\right)\mathrm{for}{u}_i\ge t_c,\end{array}& \text{(8-23)}\end{array}$

[0250] where

[0251] u_i=the actual arrival time of the train,

[0252] t_i=the time at which late penalties begin to accrue,

[0253] h_i=the size of the step penalty (in k$),

[0254] m_i=the rate of the linear portion of the penalty (in k$/hr.), and

[0255] t_c=the transition point where the cost function changes from a sigmoid to a line segment.

[0256] Shortening this cost function to the form A_i(u_i) for the nonce, and defining {right arrow over (u)}=(u₁, . . . ,u_n), we may express a cost function which accounts for the arrival times of all trains in the form $\begin{array}{cc}A\left(\stackrel{->}{u}\right)=\sum _{i=1}^{n_T}{A}_i\left(u_i\right).& \text{(8-24)}\end{array}$

[0257] But we also have the relationship $\begin{array}{cc}{u}_i=d_i+\frac{L}{s_i},& \text{(8-25)}\end{array}$

[0258] so we may consider an alternative representation of Equation 8-24 as $\begin{array}{cc}A\left(\stackrel{->}{s},\stackrel{->}{d}\right)=\sum _{i=1}^{n_T}{A}_i\left(s_i,d_i\right).& \text{(8-26)}\end{array}$

[0259] This latter form of the cost is appropriate to our search process, since that process is based on varying the components of the vector {right arrow over (s)} and {right arrow over (d)}.

[0260] Now to incorporate late arrival costs into the search, we extend the cost function of Equation 8-18 to the form

[0261] D({right arrow over (s)},{right arrow over (d)})=η₁C({right arrow over (s)},{right arrow over (d)})+η₂E({right arrow over (d)})+η₃A({right arrow over (s)},{right arrow over (d)})  (8-27),

[0262] where the η_iare weighting factors satisfying

[0263] η₁+η₂+η₃=1  (8-28).

[0264] The choices for these weights must be determined by experiment, and in one embodiment of the present invention it is possible to vary them iteratively as the search progresses. Individual users of the present invention may assign these weights as determined by the characteristics of the corridor and the costs imposed to the railroad for the various effects built into the search algorithm. In the preferred embodiment, these weights take on values as determined in conjunction with the discussion of Equation (7-8) above.

[0265] Gradient Search to Optimize Schedule Feasibility, Early Departures, and Economic Costs by Varying Only the Train Departure Times

[0266] A search using the late arrival cost function of equation 8-27 may involve variation of only the departure times {right arrow over (d)}, in which case the gradient by which the search is directed is of a form analogous to that shown in Equation 8-19. We may then express a component of the gradient vector in the form $\begin{array}{cc}\frac{\partial }{\partial d_k}\left(D\left(\stackrel{->}{s},\stackrel{->}{d}\right)\right)={\eta }_1\frac{\partial }{\partial d_k}\left(C\left(\stackrel{->}{s},\stackrel{->}{d}\right)\right)+{\eta }_2\frac{\partial }{\partial d_k}E\left(\stackrel{->}{d}\right)+{\eta }_{3}A\left(\stackrel{->}{s},\stackrel{->}{d}\right).& \text{(8-29)}\end{array}$

[0267] Borrowing from Equations 8-6 and 8-20, we expand Equation 8-29 to the form $\begin{array}{cc}\begin{array}{c}\frac{\partial }{\partial d_k}\left(D\left(\stackrel{->}{d},\stackrel{->}{s}\right)\right)={\eta }_1\left[\sum _{y_{\mathrm{ik}}\in I}\frac{\partial }{\partial y_{\mathrm{ik}}}\left(L_{\mathrm{ik}}\left(y_{\mathrm{ik}}\right)\right)\left(\frac{-s_i{s}_k}{s_i-s_k}\right)\right]-\\ {\eta }_2\left(1-\sigma \left(d_k-e_k;{\alpha }_k^{\prime }\right)\sigma \left(d_k-e_k;{\alpha }_k^{\prime }\right)\right)+\\ {\eta }_3\frac{\partial }{\partial d_k}\left(A\left(\stackrel{->}{s},\stackrel{->}{d}\right)\right)\\ ={\eta }_1\left[\sum _{y_{\mathrm{ik}}\in I}\frac{\partial }{\partial y_{\mathrm{ik}}}\left(L_{\mathrm{ik}}\left(y_{\mathrm{ik}}\right)\right)\left(\frac{-s_i{s}_k}{s_i-s_k}\right)\right]-\\ {\eta }_2\left(1-\sigma \left(d_k-e_k;{\alpha }_k^{\prime }\right)\sigma \left(d_k-e_k;{\alpha }_k^{\prime }\right)\right)+\\ {\eta }_3\frac{\partial }{\partial u_k}\left(A_k\left(u_k\right)\right)\frac{\partial }{\partial d_k}\left(u_k\right)\\ ={\eta }_1\left[\sum _{y_{\mathrm{ik}}\in I}\frac{\partial }{\partial y_{\mathrm{ik}}}\left(L_{\mathrm{ik}}\left(y_{\mathrm{ik}}\right)\right)\left(\frac{-s_i{s}_k}{s_i-s_k}\right)\right]-\\ {\eta }_2\left(1-\sigma \left(d_k-e_k;{\alpha }_k^{\prime }\right)\sigma \left(d_k-e_k;{\alpha }_k^{\prime }\right)\right)+\\ \{\begin{array}{c}{\eta }_3{h}_k{\alpha }_k\left(1-\sigma \left(u_k;{\alpha }_k\right)\right)\sigma \left(u_k;{\alpha }_k\right)\mathrm{for}{u}_k\le t_c\\ {\eta }_3{m}_k\mathrm{for}{u}_k>t_c,\end{array}\end{array}& \text{(8-30)}\end{array}$

[0268] which, with the help of Equation 8-7, provides an explicit representation of the components of the gradient of D({right arrow over (s)},{right arrow over (d)}) when only the train departure times are varied.

[0269] Gradient Search to Optimize Schedule Feasibility, Early Departures, and Economic Costs by Varying Only the Train Speeds

[0270] If the departure times of trains are held constant, and speeds are varied, then the gradient used to alter the speed vector {right arrow over (s)}=(s₁, . . . ,s_n) during the search is of the form

[0271] ∇^{({right arrow over (s)})}D({right arrow over (s)},{right arrow over (d)})=∇(η_nC({right arrow over (s)},{right arrow over (d)})+η₂E({right arrow over (d)})+η₃A({right arrow over (s)},{right arrow over (d)}))=η₁∇C({right arrow over (s)},{right arrow over (d)})+η₃∇A({right arrow over (s)},{right arrow over (d)})  (8-31),

[0272] Where E({right arrow over (d)}) is independent of the train speed. We therefore can obtain the k-th component of this gradient as $\begin{array}{cc}\begin{array}{c}\frac{\partial }{\partial s_k}\left(D\left(\overrightarrow{s},\overrightarrow{d}\right)\right)={\eta }_1\frac{\partial }{\partial s_k}\left(C\left(\overrightarrow{s},\overrightarrow{d}\right)\right)+{\eta }_3\frac{\partial }{\partial s_k}\left(A\left(\overrightarrow{s},\overrightarrow{d}\right)\right)\\ ={\eta }_1\frac{\partial }{\partial s_k}\left(C\left(\overrightarrow{s},\overrightarrow{d}\right)\right)+{\eta }_3\left[\frac{\partial }{\partial u_k}\left(A\left(\overrightarrow{u}\right)\right)\frac{\partial }{\partial s_k}\left(u_k\right)\right]\\ ={\eta }_1\frac{\partial }{\partial s_k}\left(C\left(\overrightarrow{s},\overrightarrow{d}\right)\right)+{\eta }_3·\\ \{\begin{array}{cc}\frac{-h_k{\alpha }_{k}L\left(1-\sigma \left(u_k;{\alpha }_k\right)\right)\sigma \left(u_k;{\alpha }_k\right)}{s_k^2}& \left(u_k\le t_k\right)\\ \frac{-{\mathrm{Lm}}_k}{s_k^2}& \left(u_k>t_k\right),\end{array}\end{array}& \text{(8-32)}\end{array}$

[0273] where the first term on the right side of Equation 8-32 can be expressed in a completely explicit form by reference back to Equation 8-12.

[0274] Gradient Search to Optimize Schedule Feasibility, Early Departures, and Economic Costs by Varying Both Train Departure Times and Train Speeds

[0275] In this case, both {right arrow over (d)} and {right arrow over (s)} are variables in the complete cost function of Equation 8-25, so the gradient takes the form

[0276] ∇^{({right arrow over (s)},{right arrow over (d)})}(D({right arrow over (s)},{right arrow over (d)}))=η₁∇^{({right arrow over (s)},{right arrow over (d)})}(C({right arrow over (s)},{right arrow over (d)}))+η₂∇^{({right arrow over (d)})}(E({right arrow over (d)}))+η₃∇^{({right arrow over (s)},{right arrow over (d)})}(A({right arrow over (s)},{right arrow over (d)}))  (8-33).

[0277] Again regarding the gradient in vector form, the components of the first term of the sum on the left of Equation 8-33 are readily obtainable with the help of Equations 8-14, as explicitly represented with the help of Equations 8-7, 8-11, and 8-12. The components of the second term can be obtained using Equation 8-20, and the components of the third term are obtained using Equations 8-30 and 8-32.

[0278] Including Maximum Speed Limitation Effects in the Gradient Search

[0279] A component of the cost function that would rise sharply in value as the speed s_iof a train T_ibecame close to the maximum speed s_i^(max)specified for the train was developed above. That component had the formulation (see Equation 7-9) $\begin{array}{cc}V\left(\overrightarrow{s}\right)=\sum _{i=1}^{n_T}\sigma \left(s_i-s_i^{\left(\max \right)}\right),& \text{(8-34)}\end{array}$

[0280] and occurred as a weighted term of the cost function, i.e.,

[0281] G({right arrow over (s)},{right arrow over (d)})=η₁C({right arrow over (s)},{right arrow over (d)})+η₂A({right arrow over (s)},{right arrow over (d)})+η₃E({right arrow over (d)})+η₄V({right arrow over (s)})  (8-35)

[0282] where the sum of the weights is chosen to be 1 in the preferred embodiment. Since variation of speed is independent of the departure times of trains, we have that

[0283] ∇^{({right arrow over (d)})}V({right arrow over (s)})=0  (8-36),

[0284] so constraining the search by maximum train speeds does not effect the components of the gradient obtained as partial derivatives with respect to the departure times. Relative to the gradient terms obtained as partial derivatives with respect to train speeds, we have $\begin{array}{cc}\begin{array}{c}{\nabla }^{\left(\overrightarrow{s}\right)}G\left(\overrightarrow{d},\overrightarrow{s}\right)={\eta }_1{\nabla }^{\left(\overrightarrow{s}\right)}C\left(\overrightarrow{d},\overrightarrow{s}\right)+{\eta }_2{\nabla }^{\left(\overrightarrow{s}\right)}A\left(\overrightarrow{d},\overrightarrow{s}\right)+{\eta }_4{\nabla }^{\left(\overrightarrow{s}\right)}V\left(\overrightarrow{s}\right)\\ ={\eta }_1{\nabla }^{\left(\overrightarrow{s}\right)}C\left(\overrightarrow{d},\overrightarrow{s}\right)+{\eta }_2{\nabla }^{\left(\overrightarrow{s}\right)}A\left(\overrightarrow{d},\overrightarrow{s}\right)+\\ {\eta }_4\left(\frac{\partial }{\partial s_1}\left(V\left(\overrightarrow{s}\right)\right),\dots ,\frac{\partial }{\partial s_{n_T}}\left(V\left(\overrightarrow{s}\right)\right)\right),\mathrm{and}\end{array}& \text{(8-37)}\\ \begin{array}{c}\frac{\partial }{\partial s_k}\left(V\left(\overrightarrow{s}\right)\right)=\frac{\partial }{\partial s_k}\left(\sum _{i=1}^{n_T}\sigma \left(s_i-s_i^{\left(\max \right)};\alpha ,\beta \right)\right)\\ =\frac{\alpha }{\beta }\left(1-\sigma \left(s_i-s_i^{\left(\max \right)};\alpha ,\beta \right)\right)\sigma \left(s_i-s_i^{\left(\max \right)};\alpha ,\beta \right),\end{array}& \text{(8-38)}\end{array}$