Piecewise Linear Approximation of Plane Curves with Restrictions in Computer-Aided Design of Railway Routes

Under study is approximation problems of plane curves defined by a sequence of points. These problems arise in computer-aided design planning of linear structures. The approximating curve consists of elements of a certain type (straight-line segments, arcs or 2nd degree parabolas, and also clothoids). The parameters of these elements are constrained. Moreover, the number of elements is unknown. The article deals with the simplest of the tasks of the per-element approximation, in which the elements are line segments with presence of a number of restrictions. The problem is solved by dynamic programming. An algorithm for solving the problem under the restrictions of special kind is considered. This algorithm for new CAD of railways gives decrease of construction cost when designing the optimal longitudinal profile.


Introduction
A railway line is a three-dimensional curve satisfying a number of restrictions.
The traditionally required three-dimensional curve is represented by two plane curves: horizontal and vertical alignment (later : plan and profile) ( Figure 1). The plan curve is a projection of the 3D curve to the coordinate plane XOY, and the longitudinal profile is the coordinate z as a function of the length s of the curve in the plan.
Project line is broken, which must meet a number of technical restrictions. In the current computer aided design software (CAD) the planning is carried out manually [1,2].
In the new CAD software the optimal longitudinal profile is designed by computer for each alternative in the plan [3].The project line is a broken line (polyline), which must meet a number of technical restrictions. The number of elements of the project polyline is unknown. Therefore, we believe that its nodes and the nodes of the longitudinal ground profile have the same abscissa. The longitudinal ground profile is always represented by a polyline with an uneven step between vertices, and this assumption makes it possible to fix the number of elements: n (the dimension of the problem) and s i -the length of the elements (Figure 1). Using nonlinear programming we obtain a broken line with more elements than we need, but its deviation from the required Z(s) is small because the project line is broken and must meet a number of technical restrictions [4].The idea is to find the broken line by solving the optimization problem, and then transform it into a line with elements whose length is no shorter than the predefined value, thereby defining the real dimension of the problem and the initial approach, and in the final step, perform an optimization taking into account all the restrictions and necessary revisions to the objective function. Therefore the design of the vertical alignment should be carried out in three stages. On the first and the third stages the problem is solved using a new non-linear programming algorithm [4]. The objective function is the minimum construction cost. The mathematical model takes into account the structure of the cross-section profiles of the roadbed, the presence of culverts, bridges, earth works (cuts and fills), methods of excavation, etc. In the past, a heuristic algorithm of sequential local approximation was used on the second stage. The new algorithm solves this approximation problem using dynamic programming [5,6].The presentation of this algorithm is the purpose of this article.

Statement of the Problem
Approximated curve is given by a number of points in the Cartesian coordinates L 1 = {x i , y i } (i = 1, 2, ..., n), i.e. represented as a broken line. It's required to find another broken line with the nodes L 2 = (X j , Y j ) (j = 1,2, ..., k), with the following conditions: where K (L 1 , L 2 ) predefined quantitative criterion of deviation L 2 from L 1 , and the system f r (X, Y) ≤ 0 defines the restrictions on the parameters of the elements of L 2 . The number of restrictions is defined by the number of elements (k-1), which is unknown. The first point (X 1 , Y 1 ) and the last point (X k , Y k ) are given. Additionally, the direction at start and end point may be specified (elements of contiguity). If L 1 and L 2 correspond to single -valued functions L 1 (x) and L 2 (x), then Here j is number of the nodes of the designed broken line. Additionally we can set limits at certain points in which the abscissa x r may not coincide with the abscissae of the given points. 4. L 2 (x r ) -L 1 (x r ) ≤v r and / or w r ≤L 2 (x r ) -L 1 (x r ) . Here x 1 <x r <x n . In the design of the longitudinal profile the restriction of type 1 is a discrete analogue of the restriction on the first derivative of the approximating function.
The values c j (X j ,X j+1 ) and d j (X j ,X j+1 ) depend on the j-th element position relatively to the ends of the curves in the route plan. They define the limits of the slopes.
The values p j (X j ) and q j (X j ) depend on the position of X j relatively to the transition curves (clothoids) in the route plan. They define the limits of the differences of the slopes adjacent elements of the longitudinal profile respectively over the "hump" and in the "pits". The inequality 3 together with 1 defines the minimum and maximum length of the elements. The inequality 4 limits the required extreme deviations of the sought broken line at some points.. In particular, the fixed point can be given. Values l min , l max , v r and w r are given. X j +1 -X j changed discretely, for example in steps of 10 m.
Restrictions do not allow to take the given broken line L 1 (x) as a solution.
Other criteria, in particular min max | L 2 (x)-L 1 (x)|, may be used instead of a quadratic criterion.
In the design of the longitudinal profile of railway the criterion may be the total amount of earthwork. If the soils of cuts are not used for the construction of fills, the criterion could be the cost of construction.

The Method of Solution
The unknown number of elements k does not allow to use the methods of nonlinear programming. Additionally specifying a deviation band (search area) the problem can be solved using dynamic programming.
Regarding the original curve within the search grid is broken variation (Figure 2). The transition from a given initial point A to an end point B in the final set is considered as a multi-step process. In the first step it considers all the restrictions on the allowable connections between the starting point A and the grid nodes in the range of the first element (l min ≤ X 2 -X 1 ≤ l max ), the second step -all valid connections between the end of the first element and the following nodes of the grid, so that l min ≤X 3 -X 2 ≤ l max , etc.
If for the next j-th element х n -X j+1 ≤l min , then it takes X j+1 = х n and Y j+1 = Y k . To avoid the exhaustive search, it is necessary to introduce the concept of "state of the system." In accordance with Bellman's principle [2,3] of all the acceptable paths that lead to a common state of the system, we may keep only one, for which the value of the criterion is minimal.
Because of the restrictions of type 2 we can't compare paths that converge at a single node of the grid. Therefore, as the "state of the system" for all paths we consider the last element ( but not a last node) and then we may compare the possible paths of achievement the current node for which the last element is the same. The number of elements of the paths ( trajectories) which are compared can be different. 6. Instead of element wise construction of trajectories each point inside x 1 + l min ≤ x ≤ x n -l min is considered as the beginning of the next element, no matter how many element the incoming path contains. Iteratively we construct for each point a set of admissible variants of next element ( Figure 3. "Outgoing fan").

Algorithm
At x> x n -2l min this "fan" contains only one element and its end is point B. "Outgoing fan" is built only from the points with a nonzero number of incoming elements. For each element of "outgoing fan" only one of the incoming elements is remembered. For this pair (incoming and outgoing elements) the restriction of type 2 is executed and the value of the objective function is minimal. This value is added to the estimate of the element under consideration.
Thus at each point in the x≥ x 1 + 2l min we generate a lot of incoming elements (Figure 3. "Incoming fan") and for each of them we calculate the estimate of the path which begins at point A and ends at the current point. The number of the input elements is remembered. Any "system state" is given by a pair of integers ( begin and end of the element), and for each state one integer (link), and one real number (path estimate ) are stored.
7. If the neighboring element on the right is given then at the last step when connecting to the point B we take it into account. In this case we check the restriction of type 2 in the point B and we eliminate the unsuitable elements.
8. The optimal trajectory is restored by moving in reverse.

Approximation of Multi -valued Functions
The described algorithm can be used for other purposes, including approximation of multi-valued functions.
For multi-valued functions instead of deviations on the graph ordinates we can use the deviations of the normals to the original curve. (Figure 4). At each node of the original broken line we calculate direction to the center of the circle connecting the three adjacent nodes. We identify this direction with the direction of the normal to the original line in the given node.
Calculation of the deviations of the next element of the approximating line from the original line is reduced to alternately calculating the coordinates of the intersection of this element and the normals to the original line at the nodes which are captured by this element. In Figure 4 there are the normals n1, n2, n3.
If the criterion of the quality of the approximation is the maximum deviation, then the calculation of the deviation in other points of the constructed element is not needed. Indeed, because of the linearity of elements the deviation in the intermediate points of the element from the original broken line cannot exceed the maximum deviation at the ends. Therefore, we calculate the deviations only on the normals and at the ends of the element and we take the maximum of them. (Figure 4).
Approximation of curves that are not graphs of a singlevalued functions, can be made with restrictions. Thus, the restrictions of type 4 go to the limit for deviations on the normals at some points (in particular, they can be the fixed points). Restriction of type 1 does not have much sense, and the restrictions of type 2 go to the limit on the angle at the node of the designed line. Except this we do not need to change anything in the above approximation algorithm.
In theory, as well as the approximation of single-valued functions with a slight complication of the algorithm we can solve the problem under an additional condition of adding circles with a given radius in the breaks of the designed line. Moreover, we can use several values of radius and elect the better of them. Such problems arise by projecting roads, pipelines, etc. However, in this case, we have to consider significantly larger deviations of designed line from the given line. The computational complexity is high if the number of discretes on the normals is large.

Results
It is known that the practical application of dynamic programming for a large number of states at each step involves dealing with significant computational complexity which dramatically increases with the number of states [5,6]. This phenomenon is called "curse of dimensionality." The purpose of experimental calculations was to establish the possibility of practical application of the algorithm for designing the longitudinal profile of new railways within haul on public computers with the search area width 0.5-0.6 m Due to the restrictions of type 1 and 2 we cannot perform calculations with a grid with large sells, followed by splitting them into a finer grid depending on the solution obtained. Therefore, when the width of the search area of 0.5-0.6 m we have 60 points or more at each vertical .
It was established experimentally that when the number of discretes on each vertical is from 40 to 100 and the length is from 2,500 to 3,000 there aren't any significant computational difficulties.
This enables to use the proposed algorithm not only in the second stage of the above described scheme of longitudinal profile design in conditions of a rugged relief, but also in the single stage of design under flat relief conditions. In the latter case, we design the fill whose height is not less than the given value. Therefore the grid is not divided according to the ground profile but to the appropriate line.
If we want to use dynamic programming to design of the longitudinal profile on a rugged relief without the use of nonlinear programming, then the width of the search area is equal to 10 m and more and computational difficulties are insuperable.

Discussion
The method of solving the problem is based on the pioneering work of R. Bellman, who 50 years ago solved the problem of piecewise linear approximation of planar curves [7]. We used the key concepts of "step process", the "system state", etc. as is common in his works [5][6][7]. The difference is in presence of restrictions (1-4 types). These differences demanded to change the formalization of these concepts ( "step process" and "state of the system") and led to a significant complication of the problem.
The proposed algorithm is one element of a multi-stage technology computer-aided design of the longitudinal profile of the railways. This technology is based on complex optimization methods and has clear advantages over the interactive elaboration of design solutions.
The use of complex mathematical models can solve design problems with their relationship and create intelligent CAD. Designs of longitudinal profile using optimization algorithms give objectively comparison the various route options and choose the most appropriate one. This approach is more promising than the various heuristic algorithms for railways routes design.
Practical tasks were solved on old computer with RAM =512 Mb, CPU 2.0 GHz . For length of the longitudinal profile 25 km the process time was a few minutes. Therefore practical tasks may be solved for a suitable time on modern computer.
The algorithm approximation of discretely defined plane curves with restrictions may find other applications that have nothing to do with designing of railways routes .

Conclusion
To apply dynamic programming it is necessary to identify the key concepts: the "step process", the "system state" and give the rule how to construct a set of states at each step, how to reject unpromising paths leading to the current state, and what must be stored to recover the optimal trajectory.
All this is done in this article with respect to the piecewise linear approximation of the plane curve given sequence of points.
Due to restrictions we can not take a single point on the vertical as the "system state". It is necessary to add the value of slope (derivative) in it. Search of the optimal number of elements and their lengths with all the restrictions led to a significant complication of the problem. But such search provides improvement of the quality of design solutions, because any small change in the alignment, will result in a significant change in the total cost [8,9]. This is the main difference a new CAD software from previously developed [10][11][12].
The described algorithm is used in two CAD systems. The first is for the design of new railway lines, and the secondfor the design of the reconstruction of railways. A new algorithm assumes the use of a piecewise parabolic approximation in CAD of roads. This algorithm will be published in a future article.