=Paper=
{{Paper
|id=Vol-1987/paper19
|storemode=property
|title=Decommissioning of Nuclear Facilities: Minimum Accumulated Radiation Dose Routing Problem
|pdfUrl=https://ceur-ws.org/Vol-1987/paper19.pdf
|volume=Vol-1987
|authors=Alexander G. Chentsov,Alexey M. Grigoryev,Alexey A. Chentsov
}}
==Decommissioning of Nuclear Facilities: Minimum Accumulated Radiation Dose Routing Problem==
https://ceur-ws.org/Vol-1987/paper19.pdf
Decommissioning of Nuclear Facilities: Minimum
Accumulated Radiation Dose Routing Problem
Alexander G. Chentsov Alexey M. Grigoryev Alexey A. Chentsov
Krasovskii Institute of Krasovskii Institute of Krasovskii Institute of
Mathematics and Mathematics and Mathematics and
Mechanics UB RAS Mechanics UB RAS Mechanics UB RAS
S.Kovalevskaya Str. 16, S.Kovalevskaya Str. 16, S.Kovalevskaya Str. 16,
620990 Yekaterinburg, 620990 Yekaterinburg, 620990 Yekaterinburg,
Russia. Russia. Russia.
chentsov@imm.uran.ru ag@uran.ru chentsov a a@mail.ru
Abstract
The problem that we consider has its prototype in an engineering prob-
lem that arises in dismantling, one at a time, a system of elements emit-
ting hazardous radiation. Its possible applications can be connected
with optimizing the actions of nuclear power plant staff in emergen-
cies or disaster cleanup operations (e.g. those comparable with Cher-
nobyl or Fukushima Daiichi). The main features of the arising routing
problem as compared with the well-known Traveling salesman problem
(TSP) are connected with the dependence of travel and dismantlement
cost functions on the set of pending tasks (not yet dismantled radi-
ation sources) and precedence constraints. For example, a radiation
source’s influence on the agent is computed by integrating the depen-
dence, which is inversely proportional to the square of the distance;
only the sources that are not yet dismantled at the time of a movement
are considered. Thus, the agent’s exposure to radiation or other haz-
ards depends on his specific trajectory and the order of dismantlement
operations.
1 Introduction
This paper is devoted to a study of an important class of discrete optimization problems, namely, the problems
concerning multiple movements under constraints. The prototype of the considered problem is the well-known
intractable traveling salesman problem (TSP), see [Garey, 1979, Ch.3]; to name just a few works on TSP,
[Gutin, 2002, Cook, 2012, Melamed, 1989, Litl, 1965, Bellman, 1962, Held, 1962]. However, in applications, one
often has to satisfy certain additional constraints. In particular, constraints arise in the problem of minimizing
the exposure of nuclear power plant staff conducting operations related to dismantling radiation sources. A
typical characteristic of the latter problem is dependence on the set of pending tasks: only the sources that are
not yet dismantled at the time of each operations do “radiate.” There are also other applications [Leon, 1996,
Copyright ⃝
c by the paper’s authors. Copying permitted for private and academic purposes.
In: Yu. G. Evtushenko, M. Yu. Khachay, O. V. Khamisov, Yu. A. Kochetov, V.U. Malkova, M.A. Posypkin (eds.): Proceedings of
the OPTIMA-2017 Conference, Petrovac, Montenegro, 02-Oct-2017, published at http://ceur-ws.org
123
Alkaya, 2010, Kinable, 2016]. Another peculiarity, which arises, for example, in the known problem of dismantling
a decommissioned nuclear power generation unit, consists in precedence constraints on the sequence of operations;
this is formalized by specifying ordered pairs (OP) of operations where the first component of the pair must
be conducted before the second one. Finally, in contrast with “ordinary” TSP, in applications, the “cities”
visited may possess some internal structure and are thus rendered not as “cities” but as clusters of cities, or
megalopolises, which introduces a certain additional variation into the movements—the problem develops a two-
layer hierarchy: layer one is the sequencing of megalopolises through a permutation of their indices—the route,
and layer two is the choice of a “track” along the given route, i.e., the exact versions of movements through
the megalopolises. The arguments leading to a two-layer optimization problem are stated in [Chentsov, 2008],
and the applications of methods developed in [Chentsov, 2008] for problems of nuclear power generation are
discussed in [Korobkin, 2012, Tashlykov, 2011]. In this paper, we use a rather general dynamic programming
(DP) procedure.
2 Problem Statement and Discussion
△
Let us start with general definitions, notions, and notations. The symbol = denotes equality by definition. A
set, all elements of which are themselves sets is called a family. The real line is denoted by R; in addition,
△ △ △
R+ = {ξ ∈ R | 0 6 ξ}, N = {1; 2; . . .} and N0 = N ∪ {0} = {0; 1; 2; . . .}, N ⊂ N0 ⊂ R. For p ∈ N0 and q ∈ N0 ,
△
assume p, q = {t ∈ N0 | (p 6 t) & (t 6 q)}. By R+ [S] we denote the set of all functions from a nonempty set S
into R+ , i.e., the set of all nonnegative real functions on S. For every ordered pair (OP) z = (a, b) of arbitrary
objects a and b, denote by pr1 (z) and pr2 (z), respectively, its first and second elements: pr1 (z) = a, pr2 (z) = b.
For an object x, the corresponding singleton set is denoted by {x}. For arbitrary objects a, b, and c, by
△
convention [Dieudonne, 1960, p.17], assume (a, b, c) = ((a, b), c); thus, a triple is a certain ordered pair. For three
△
arbitrary sets A, B, and C, following [Dieudonne, 1960, p.17], we assume A × B × C = (A × B) × C, thus, for
µ ∈ A × B and ν ∈ C, we have (µ, ν) ∈ A × B × C. For a set H, denote by P(H) the set of all subsets of H,
′ △
i.e., the power set of H; P (H) = P(H) \ {∅}; denote by Fin(H) the family of all (nonempty) finite sets of
′ ′
P (H). If H is a finite set, then Fin(H) = P (H). To a nonempty finite set K, assign its cardinality |K| ∈ N
(the number of its elements) and the (nonempty) set (bi)[K] of all bijections [Cormen, 1990, p.87] of the set
△
1, |K| = {j ∈ N | j 6 |K|} onto K; in addition, |∅| = 0. A permutation of a nonempty set T is [Cormen, 1990,
p.87] a bijection of T onto itself; for every permutation α of T , there exists an inverse permutation α−1 of the
set T : α(α−1 (t)) = α−1 (α(t)) = t ∀t ∈ T . The permutations of indices used below are interpreted as sequences
of visiting the objective sets—the megalopolises.
Fix a nonempty set X, a point x◦ ∈ X, a number N ∈ N, N > 2, (nonempty finite) sets M1 ∈
′ ′
Fin(X), . . . , MN ∈ Fin(X), and relations M1 ∈ P (M1 × M1 ), . . . , MN ∈ P (MN × MN ). Here and below,
assume (x◦ ∈ / Mj ∀j ∈ 1, N ) & (Mp ∩ Mq = ∅ ∀p ∈ 1, N ∀q ∈ 1, N \ {p}). The objects M1 , . . . , MN are
regarded as the megalopolises to be visited. For each j ∈ 1, N , OPs from Mj determine the feasible interior jobs
connected with visiting Mj : the first element of such an OP is the entry point, and the second element is the
exit point. We consider the processes of the following [Chentsov, 2013, (3.3)] form:
x◦ → (pr1 (z1 ) ∈ Mα(1) pr2 (z1 ) ∈ Mα(1) ) → . . . → (pr1 (zN ) ∈ Mα(N ) pr2 (zN ) ∈ Mα(N ) ), (2.1)
where z1 ∈ Mα(1) , . . . , zN ∈ Mα(N ) , and α is a permutation of indices of 1, N , hereinafter referred to as route.
Assume that in (2.1) the straight arrows denote the exterior movements and the wavy arrows denote the motions
△
connected with (interior) jobs. The objects of our choice
△
(∪ are α,)z1 , . . . , zN . Let△ Mj = {pr z ∈ Mj } ∈
(∪2 (z) : )
◦ N ◦ N
Fin(Mj ) ∀j ∈ 1, N . In view of this, assume X = {x } ∪ i=1 Mi ∈ Fin(X), X = {x } ∪ i=1 Mi ∈ Fin(X).
In connection with (2.1), note that both exterior movements and interior jobs conducted as the megalopolises
are visited are measured through the given functions, and the results of these measurements are aggregated
△
additively, which is the natural way for many applications. Here and below, assume P = (bi)[1, N ], and also
assume that each specific α ∈ P (see (2.1)) must satisfy the precedence constraints, which are defined by means
of the set K ∈ P(1, N × 1, N ) ; OPs from K will be called address pairs; the feasibility of α ∈ P can be reduced
to the following requirement: for z = (i, j) ∈ K, the set Mi must be visited before Mj . Here and below, assume
′
[Chentsov, 2013, (3.11)] the following: ∀K0 ∈ P (K)∃z0 ∈ K0 : pr1 (z0 ) ̸= pr2 (z) ∀z ∈ K0 . Then, the set A of
124
all feasible (with respect to precedence) routes from P has the following form [Chentsov, 2013, (3.12)]:
△ ′
A = {α ∈ P|α−1 (pr1 (z)) < α−1 (pr2 (z)) ∀z ∈ K} ∈ P (P). (2.2)
Note that pr1 (z) ̸= pr2 (z) ∀z ∈ K. As evident from (2.1), the choice of some route α ∈ A does not yet
determine the whole process, it must be supplemented with some track z1 , . . . , zN . A track has to agree
with the given route, and the solution is constructed in the form of OP consisting of a route and a track,
where the route must be an element of A. It is useful to bear in mind that such interpretation allows one to
consider not only (2.1) but also the partial processes connected with visiting megalopolises Mk , k ∈ K, where
K ⊂ 1, N . Such interpretation is also useful in connection with the dynamic programming (DP); it was laid
out in [Chentsov, 2013, Chentsov, 2015, Chentsov, 2016]. In this paper, we omit its details and devote most of
our attention to the computational concerns. Thus, we limit ourselves to a description of “complete” routes and
“tracks” (trajectories), similar to (2.1). Thus, if α ∈ P, denote by Zα the set of all tuples (zi )i∈0,N : 0, N → X×X
such that (z0 = (x0 , x0 )) & ( zt ∈ Mα(t) ∀t ∈ 1, N ); Zα is a nonempty finite set. Feasible solutions will be
△ ′
represented by the OPs (α, (zi )i∈0,N ), where α ∈ A and (zi )i∈0,N ∈ Zα . Let N = P (1, N ).
Consider the following cost functions (they are assumed to be known): c ∈ R+ [X × X × N], c1 ∈ R+ [X × X ×
N], . . . , cN ∈ R+ [X × X × N], f ∈ R+ [X]. In terms of these functions, let us define the additive criterion: for
α ∈ A and (zi )i∈0,N ∈ Zα , set
△ ∑
N
Cα [(zi )i∈0,N ] = [c(pr2 (zt−1 ), pr1 (zt ), {α(j) : j ∈ t, N })+
t=1 (2.3)
+ cα(t) (zt , {α(j) : j ∈ t, N })] + f (pr2 (zN )).
Our principal problem is as follows:
Cα [(zi )i∈0,N ] → min, α ∈ A, (zi )i∈0,N ∈ Zα ; (2.4)
for each problem (2.4), there exists a nonempty set of optimal feasible solutions and the problem’s value—the
following extremum:
△
V = min min Cα [(zi )i∈0,N ] ∈ R+ . (2.5)
α∈A (zi )i∈0,N ∈ Zα
Our goal is to find the (global) extremum V (2.5) and some optimal feasible solution (α0 , (zi0 )i∈0,N ), where
α ∈ A and (zi0 )i∈0,N ∈ Zα0 ; evidently, they satisfy Cα0 [(zi0 )i∈0,N ] = V .
0
3 Dynamic Programming. Layers of Bellman Function
Problem (2.4) is solved by means of a variety of DP [Chentsov, 2013, Chentsov, 2015, Chentsov, 2016], which
we specify in a brief form. We use the construction from [Chentsov, 2013, Chentsov, 2015, Chentsov, 2016],
which is based on the layer structure of the Bellman function; the algorithm of their construction is sketched
below. Let us consider the construction of state space layers by means of the crossing-out operator (for tasks
△
in the task list) [Chentsov, 2008, Pt. 2]: I : N → N, specifically, for K ∈ N, set Ξ(K) = {z ∈ K | (pr1 (z) ∈
△
K) & (pr2 (z) ∈ K)} and I(K) = K \ {pr2 (z) : z ∈ Ξ(K)}. This operator is used to construct the layers of
△
the state space. To this end, consider a set G = {K ∈ N | ∀z ∈ K (pr1 (z) ∈ K) ⇒ (pr2 (z) ∈ K)} clearly,
△
1, N ∈ G. We sort the mentioned lists by their cardinality, Gs = {K ∈ G | s = |K|} ∀s ∈ 1, N . Then, the family
△ △
{Gj : j ∈ 1, N } is a partition of G. In addition, G1 = {{t} : t ∈ 1, N \ K1 }, where K1 = {pr1 (z) : z ∈ K}. Also
[Chentsov, 2013, Chentsov, 2015, Chentsov, 2016], Gs−1 = {K \ {j} : K ∈ Gs , j ∈ I(K)} ∀s ∈ 2, N . Thus we
obtain a recurrence procedure GN → GN −1 → . . . → G1 . Based on this procedure, we construct the state space
layers, which we denote D0 , D1 , . . . , DN . Specifically, D0 = {(x, ∅) : x ∈ M}, where M is by definition the
△
union of all the sets Mi , i ∈ 1, N \ K1 . Next, DN = {(x0 , 1, N )} (the singleton that contains the OP (x0 , 1, N )).
If s ∈ 1, N − 1 and K ∈ Gs , we have [Chentsov, 2015] the following sequence of definitions:
△ △ ∪
Js (K) = {t ∈ 1, N \ K | {t} ∪ K ∈ Gs+1 }, Ms [K] = Mj ,
j∈Js (K)
△ ′
Ds [K] = {(x, K) : x ∈ Ms [K]} ∈ P (X × Gs ).
125
Thus, for s ∈ 1, N − 1, the layer Ds is defined as the union of all sets Ds [K], K ∈ Gs . It is easy to see that
D0 ̸= ∅, D1 ̸= ∅, . . . , DN ̸= ∅. If s ∈ 1, N , (x, K) ∈ Ds , j ∈ I(K), and z ∈ Mj , then (pr2 (z), K \ {j}) ∈ Ds−1 .
Recurrence procedure for construction of the layers. Consider a system of funcions v0 ∈ R+ [D0 ], v1 ∈
△
R+ [D1 ], . . . , vN ∈ R+ [DN ]. First, define v0 ∈ R+ [D0 ] by the condition v0 (x, ∅) = f (x) ∀x ∈ M. Further
constructions implement the following recurrence scheme: for s ∈ 1, N , if the function vs−1 ∈ R+ [Ds−1 ] is
already constructed, vs ∈ R+ [Ds ] can be determined by the rule
△
vs (x, K) = min min [c(x, pr1 (z), K) + cj (z, K) + vs−1 (pr2 (z), K \ {j})| ∀(x, K) ∈ Ds . (3.1)
j∈I(K) z∈Mj
From the general constructions of [Chentsov, 2013, Chentsov, 2015, Chentsov, 2016], there follows the equality
V = vN (x0 , 1, N ), (3.2)
The procedure v0 → v1 → . . . → V, which is completely defined by virtue (3.1), can be regarded as an algorithm
for determining V , during which, only one layer of the Bellman function is retained in the computer’s RAM; for
details, refer to [Chentsov, 2016.1]. Note that a similar idea was also proposed in [Lawler, 1979]. Construction
of optimal solution, corresponds to [Chentsov, 2013, Chentsov, 2016].
4 Scheme of Independent Computations for Layers of Bellman Function
Further exposition follows a version of the general construction [Chentsov, 2013, Chentsov, 2012,
Chentsov, 2012.1] connected with conducting independent computations. Recall that (see (3.2)) states the “fi-
nal” formula for determining V through vN −1 . Thus, further constructions require one to determine vN −1 ,
which is conducted by the procedure below. The function vN −1 is defined on the set DN −1 , which is de-
fined through GN −1 . Here and below, assume N > 3; we have GN −1 = {1, N \ {j} : j ∈ I(1, N )} since
GN = {1, N }, it is possible to determine DN −1 . Assume all the layer-functions will be constructed on n
△
nodes, where (here and below) n = |GN −1 | = |I(1, N )| ∈ N. To facilitate these constructions, each layer
D0 , D1 , . . . , DN will be interpreted as a union of n subsets; each such subset will be distributed to a specific
node, which would construct its fragment of the layer. This construction will be conducted through auxiliary
discrete dynamic systems (DDS), the trajectories of which will be determined through systems of inclusions.
Thus, assume K ∈ GN −1 , and let T[K] be the set of all the trajectories (Kt )t∈0,N −2 : 0, N − 2 → G such that
△
(K0 = K) & (∀τ ∈ 1, N − 2 ∃s ∈ I(Kτ −1 ) : Kτ = Kτ −1 \ {s}). Trajectories of the system are not uniquely
e
determined for each given K; thus, for each “time” t ∈ 0, N − 2, we have a reachability set (RS) T[K; t], defined
to be the set of all the lists Kt such that (Kτ )τ ∈0,N −2 are the trajectories from T[K]. Now, the families of
feasible lists of fixed cardinality can be determined through these RSs, specifically,
∪
GN −(t+1) = e
T[K; t] ∀t ∈ 0, N − 2. (4.1)
K∈GN −1
Let us also note that the mentioned RSs could [Chentsov, 2012.1, Proposition 16] be determined in a recurrent
fashion for a fixed condition: if K ∈ GN −1 and t ∈ 0, N − 3, then
e
T[K; e
t + 1] = {P \ {h} : P ∈ T[K; t], h ∈ I(P )}. (4.2)
Through RSs, we construct the “individual” (to each computation node—its own) state space layers: for K ∈
GN −1 and s ∈ 1, N − 1, set
△ ∪ ′
Ds [K] = Ds [P ] ∈ P (Ds ). (4.3)
e
P ∈T[K;N −(s+1)]
These layers (4.3) possess the property,
(pr2 (z), Q \ {s}) ∈ Dl [K] ∀(x, Q) ∈ Dl+1 [K] ∀s ∈ I(Q) ∀z ∈ Ms . (4.4)
In view of (4.3), for K ∈ GN −1 and s ∈ 1, N − 1, assume
△
Ws [K] = (vs (x, P ))(x,P )∈Ds [K] = (v(x, P ))(x,P )∈Ds [K] ∈ R+ [Ds [K]]. (4.5)
126
△
Suppose that M0 [K] = {h ∈ 1, N \ K1 |{h} ∈ T[K;e N − 2]}. Next, for K ∈ GN −1 and h ∈ M0 [K], we have
e
{h} ∈ T[K; N − 2], where h ∈ 1, N \ K1 ; thus, {h} ∈ G1 and J1 ({h}) = {t ∈ 1, N \ {h}|{t; h} ∈ G2 }, then
∪
M1 [{h}] = Mj . (4.6)
j∈J1 ({h})
For K ∈ GN −1 , h ∈ M0 [K], and x ∈ M1 [{h}], we have the following property: (x, {h}) ∈ D1 [K] and
W1 [K](x, {h}) = min [c(x, pr1 (z), {h}) + ch (z, {h}) + f (pr2 (z))]. (4.7)
z∈Mh
So, through (4.7), we can fully determine the function W1 [K], while, to use (4.7) we have to, knowing
K ∈ GN −1 and h ∈ M0 [K], construct J1 ({h}) and then determine M1 [{h}] by means of (4.6), after which we
can finally use formula (4.7) for x ∈ M1 [{h}]. Thus we determine all the functions W1 [K], K ∈ GN −1 . Further
constructions will concern one single computation node.
A transformation of a function Wl [K] into Wl+1 [K], where K ∈ GN −1 and l ∈ 1, N − 1, is determined by a
relation similar to (3.1), which makes use of (4.4), (indeed, by (4.4) and (4.5), for (x, Q) ∈ Dl+1 [K], s ∈ I(Q),
and z ∈ Ms , the value Wl [K](pr2 (z), Q \ {s}) ∈ R+ is defined and can be used to compute Wl+1 [K](x, Q)). Thus
(see (3.1) and (4.5)), for K ∈ GN −1 , l ∈ 1, N − 2, and (x, Q) ∈ Dl+1 [K], we have [Chentsov, 2013, (10.17)]
Wl+1 [K](x, Q) = min min [c(x, pr1 (z), Q) + cj (z, Q) + Wl [K](pr2 (z), Q \ {s})]. (4.8)
s∈I(Q) z∈Mj
Fix K ∈ GN −1 . For this set, we have the corresponding layers Ds [K] (4.3), s ∈ 1, N − 1, which form a
nonempty set in the state space. For s ∈ 1, N − 1, the function Ws [K] ∈ R+ [Ds [K]] is defined. In particular, the
function W1 [K] ∈ R+ [D1 [K]] is defined; to determine its values, (4.7) should be used in view of the representation
of D1 [K]. Further construction of functions W1 [K], . . . , WN −1 [K] is conducted based on a recurrence procedure
based on (4.8). Specifically, for l ∈ 1, N − 2, the transformation of Wl [K] into Wl+1 [K] is determined by (4.8)
under K = K. Since the choice of l was arbitrary, we obtain the recurrence procedure W1 [K] → W2 [K] → . . . →
WN −1 [K]. All these computations are conducted by a single computation node independently of other nodes.
5 Construction of Layers of Bellman Function
Note the positions of [Chentsov, 2012.1, Section 7]. In particular, from [Chentsov, 2012.1, Proposition 17]: we
have the equalities ∪
Ds = Ds [K] ∀s ∈ 1, N − 1. (5.1)
K∈GN −1
Equalities (5.1) provide for a “union” of the processes conducted by the separate computation nodes. Our
intention in this is to construct the layers v1 , . . . , vN −1 . For every s ∈ 1, N − 1, we have the functions Ws [K], K ∈
GN −1 , the domains of which (that is, the sets Ds [K], K ∈ GN −1 ) form, in view of (5.1), a cover of Ds . The
above-mentioned functions agree in view of (4.5): for K1 ∈ GN −1 , K2 ∈ GN −1 , (x, P ) ∈ Ds [K1 ] ∩ Ds [K2 ],
Ws [K1 ](x, P ) = vs (x, P ) = Ws [K2 ](x, P ). (5.2)
In view of (5.1) and (5.2), we obtain the following simple rule of construction of the function vs for s ∈ 1, N − 1.
Having “particular” functions, to determine vs (x0 , K0 ), where (x0 , K0 ) ∈ Ds , first, in view of (5.1), obtain the
set K0 ∈ GN −1 such that (x0 , K0 ) ∈ Ds [K0 ]. Then, in view of (4.5), vs (x0 , K0 ) = Ws [K0 ](x0 , K0 ). Property (5.2)
means that a specific choice of K0 ∈ GN −1 with the property (x0 , K0 ) ∈ Ds [K0 ] can be arbitrary. So, all values of
the function vs can be determined and, therefore, the function itself. Thus, being in possession of functions (4.5),
we obtain all the layers v0 , v1 , . . . , vN of the Bellman function. Based on these layers, through procedures that
solve the local problems, we determine the optimal feasible solution. To construct such a solution, the computer
should retain all the layers v1 , . . . , vN of the Bellman function in its memory.
Algorithm for determining the global extremum. Consider a procedure for determining V (3.2), during
the implementation of which, in the memory of each computational node, there is only retained one “particular”
function of the form (4.5), i.e., an “individual” Bellman function layer.
Without loss of generality, consider the computational node that processes the set K ∈ GN −1 . The ultimate
goal of the procedure running on this node is to construct the function WN −1 [K] ∈ R+ [DN −1 [K]].
127
Principal steps of the iteration procedure.
1)Determine W1 [K] ∈ R+ [D1 [K]] from (4.7).
2)For s ∈ 1, N − 2, assume the “particular” function Ws [K] ∈ R+ [Ds [K]] is already constructed. Then,
through (4.8), compute the values of the function Ws+1 [K], which only use the values of the function Ws [K]
(in (4.8), set l = s). This yields the “particular” function Ws+1 [K] ∈ R+ [Ds+1 [K]].
Next, the memory holding the values Ws [K] is released, and then filled with the values of Ws+1 [K]: the
“particular” functions are overwritten.
3)After consecutively conducting step 2), we obtain WN −1 [K].
After steps 1)–3) are completed at each computation node, there are obtained all the functions
WN −1 [K], K ∈ GN −1 . (5.3)
∪
Now, construct vN −1 through the equality DN −1 = DN −1 [K]. To this end, for every state (x̂, P̂ ) ∈
K∈GN −1
DN −1 , determine, the set K̂ ∈ GN −1 such that (x̂, P̂ ) ∈ DN −1 [K̂]. Then (see (5.3)), we have the value
WN −1 [K̂](x̂, P̂ ) ∈ R+ for which, in view of (4.5), there holds the equality vN −1 (x̂, P̂ ) = WN −1 [K̂](x̂, P̂ ). Thus,
we obtain all the values of the function vN −1 .
Through the function vN −1 (the layer of the Bellman function), we determine the global extremum V . In
connection with this procedure, note [Chentsov, 2016.1], where it was proved that one can determine the value
of V without constructing the optimal feasible solution for problems of greater dimension.
6 Model and Computational Experiment.
Let us consider a specific version of the general problem; assume X = R × R (the problem on a plane). We
study a formulation where each megalopolis is a system of entries and exits for a certain zone of heightened
intensity of certain harmful factors (in particular, radiation). The necessary activities are conducted by a point
agent in the zone. The aim of activities inside a megalopolis is to go from the entry point to the radiation
source, dismantle it, and then leave the megalopolis through an exit point. The radiation left after dismantling
a source is assumed to be negligible. Thus, a megalopolis is, essentially, a “near zone” of the radiation source.
Other (non-dismantled) radiation sources are also assumed to have a nonnegligible effect both during movement
between megalopolises and during the near zone activities (interior jobs). The effect of radation is assumed to
be cumulative, whence the additive cost aggregation assumed in this paper. Another feature of the problem
statement is the fact that when moving to a source to dismantle it, as the agent enters the near zone, the agent is
assumed to be affected by it; however, after dismantlement, as the agent moves to exit the zone, the agent is not
affected by the (dismantled) source anymore. Consequently, the degree of radiation exposure of the agent during
the first stage depends on the length of the stage and, therefore, on the distance between the entry point and
the radiation source. During exterior movements, the effect of each single source is not as pronounced, however,
the cumulative effect of non-dismantled sources can not be ignored.
Passing by a radiation source s and dismantling it. Let us connect the system of megalopolises with
a system of point “radiating” objects (zi )i∈1,N : 1, N → X with the following property: zi ∈ / Mi ∀i ∈ 1, N .
Dismantlement of the objects z1 , . . . , zN is the aim of visiting the megalopolises: specifically, during a visit
to Mα(i) , the agent has to enter at some entry point pr1 (zi ), move to the point zα(i) , dismantle the radiation
source number α(i), and then reach an exit point pr2 (zi ). This implements scheme (2.1).
Let us sketch a description of the effect of a single source during the movement from the given point of the
plane to another point (as mentioned before, the exposure suffered during these movements is summed up) in
the “regular” case, when the trajectory does not pass through any “active” radiation sources. Thus, as the agent
moves from a point i to a point j, the losses (in our model), or rather, the exposure suffered during the movement
due to the effect of the (active) source s is as follows:
∫T ∫ρi,j
γs γs dρ
ci,j [{s}] = dt = = 2ρi,j ( )2 ( ), (6.1)
ρ2s,t (t) v 2ρi,j ρ + ρ2j,s − ρ2i,s − ρ2i,j + 4ρ2i,j ρ2i,s − (ρ2j,s − ρ2i,s − ρ2i,j )2
0 0
where ρ denotes the Euclidean distance (where necessary, the mentioned distances are subscripted); γs is the
intensity of the source s, and v is the movement speed. To calculate the integral use one of the following table
formulas: ∫ ∫
dR 1 R dR 1 A+R
= arctg + C, = ln +C,
2
A +R 2 A A A −R
2 2 2A A − R
128
depending on whether the sign is 4ρ2i,j ρ2i,s − (ρ2j,s − ρ2i,s − ρ2i,j )2 . If the source s lies on the trajectory of the motion
from point i to point j, the cost ci,j [{s}] is assumed to be a very large number (roughly speaking, ci,j [{s}] = ∞;
in the actual computation, the number that is several times greater than the most “costly” motion is sufficient).
Approaching radiation source s for dismantlement (determining the cost of interior jobs). In this
case, assume that the point j coincides with s and consider the exposure suffered by the agent during the approach
to the source that is to be dismantled. Thus, we still follow the model where the exposure is inversely proportional
to the square of the distance to the source; however, to avoid the ill-posedness as the agent reaches the source,
which could lead to a division by 0, the denominator of the integrand in (6.1) is increased by 1. To account for a
more intense radiation in the near zone, we add the factor 3 (which describes the case of radiation intensity not
being weakened by obstacles, the latter being possible during exterior movements). So, in this case, the losses
( ) ρ∫i,j dρ ( γs )
(the exposure) are calculated by the following formula: ci,j [{s}] = ci,s [{s}] = 3 γvs ρ2 +1 = 3 v arctg(ρi,j ).
0
Computational experiment. Let us consider the model instances of the routing problem for dismantling
radiation sources on a plane. The megalopolises, which imitate the possible entry and exit points of the chambers
housing the radiation sources, are obtained by discretizing circles: at each circle, there are 12 points, equally
spaced starting with the point with the angle coordinate of 0. Each megalopolis is assigned a point object, which
imitates the radiation source of the chamber. Let the starting point of the dismantlement process coincide with
the point of origin, i.e., x0 = (0, 0); after dismantling all radiation sources, the agent must return to the depot
(the point of origin). Recall that the function ρ is essentially the Euclidean distance. Assume the speed of the
agent is 4 times greater outside the chambers than it is within—this models the intrinsic difficulty of moving
inside the megalopolises, which is due to the presence of hardware, various structures, or mechanisms, which
hamper the rapid movement inside.
Here are some results concerning the calculations on supercomputers URAN. For the model with 30 mega-
lopolises and 30 address pairs (51 pairs in transitive closure [Schmidt, 1993]) that determine the precedence
constraints, the following results were obtained: total exposure (total losses): 222.9, total computation time:
17m. 46s. For the model with 31 megalopolis and 34 address pairs (63 pairs in transitive closure[Schmidt, 1993]),
the following results were obtained: total exposure (total losses): 226.5, total computation time: 15m. 56s. Let
us now consider the results of the accounts of the same tasks on a PC. For the model with 30 megalopolises
and 30 address pairs, the following results were obtained: total exposure (total losses): 222.9, total computation
time: 7h. 26m. 7s. For the model with 31 megalopolis and 34 address pairs, the following results were obtained:
total exposure (total losses): 226.5, total computation time: 6h. 29m. 55s. A decrease in computation time for
the problem with the greater number of megalopolises is due to the more “strict” precedence constraints, which
significantly decrease both the volume of computations involved in obtaining the values of the Bellman function
and the memory usage. For more information on the influence of precedence constraints on computation time
and memory usage of dynamic programming solutions, refer to [Steiner, 1990].
Acknowledgements
This work was supported by Russian Science Foundation (project no. 14-11-00109).
References
[Alkaya, 2010] Alkaya, A.F., Duman, E. (2010). A new generalization of the traveling salesman problem. Appl.
Comput. Math. 9(2), 162–175.
[Bellman, 1962] Bellman, R. (1962). Dynamic programming treatment of the travelling salesman problem. J.
Assoc. Comput. Mach. 9, 61–63.
[Chentsov, 2008] Chentsov, A.G. (2008). Extreme Problems of Routing and Tasks Distribution: Regular and
Chaotic Dynamics. Izhevsk Institute of Computer Research, 240 p. (in Russian).
[Chentsov, 2012] Chentsov, A.G. (2012). One parallel procedure for the construction of the Bellman function in
the generalized problem of the courier with the inner workings. Autom. Remote Control 3, 134–149.
[Chentsov, 2012.1] Chentsov, A.G. (2012). One parallel procedure for the construction of the Bellman function
in the generalized problem of the courier with the inner workings. Vestnik YuUrGU. Gos. Univ. Ser.
Mat. Model Program. 12, 53–76.
129
[Chentsov, 2013] Chentsov, A.G. (2013). To question of routing of works complexes. Vestnik Udmurtskogo Uni-
versiteta. Matematika. Mekhanika. Kompyuternye Nauki, (1), 59–82.
[Chentsov, 2015] Chentsov, A.G., Chentsov, A.A. (2015). Route problem with constraints depending on a list of
tasks. Doklady Mathematics, 92(3), 685–688.
[Chentsov, 2016] Chentsov, A.G., Chentsov, P.A. (2016). Routing under constraints: problem of visit to mega-
lopolises. Autom. Remote Control, 77(11), 1957–1974.
[Chentsov, 2016.1] Chentsov, A.G., Chentsov, A.A. (2016). On the problem of obtaining the value of routing
problem with constraints. J. Autom. Inf. Sci. 6, 41–54.
[Chentsov, 2016.2] Chentsov, A.G., Grigoryev A.M. (2016). Dynamic programming method in the route problem:
the scheme of independent calculations. Mekhatronika, Avtomatizatsiya, Upravlenie, 17(12), 834–846.
[Cook, 2012] Cook, W.J. (2012). In Pursuit of the Traveling Salesman. Mathematics at the limits of computation.
Princeton University Press, New Jersey, p. 248.
[Cormen, 1990] Cormen, T.H., Leizerson, C.E., and Rivest, R.L. (1990). Introduction to Algorithms. Cambridge:
MIT Press.
[Dieudonne, 1960] Dieudonne, J. (1960). Foundations of Modern Analysis. New York: Academic.
[Garey, 1979] Garey, M.R., Johnson, D.S. (1979).Computers and Intractability: A Guide to the Theory of NP-
Completeness. N.Y., W.H. Freeman & , p. 416
[Gutin, 2002] Gutin, G., Punnen, A.P. (2002). The Traveling Salesman Problem and Its Variations. Springer,
New York.
[Held, 1962] Held, M., Karp, R.M. (1962). A dynamic programming approach to sequencing problems. J. Soc.
Ind. Appl. Math. 10(1), 196–210.
[Kinable, 2016] Kinable J., Cire A., van Hoeve W. J. (2016). Hybrid optimization methods for time-dependent
sequencing problems. European Journal of Operational Research, Available online 24 November 2016,
ISSN 0377-2217
[Korobkin, 2012] Korobkin, V.V., Sesekin, A.N., Tashlykov, O.L., Chentsov, A.G. (2012). Routing Methods and
Their Applications to the Enhancement of Safety and Efficiency of Nuclear Plant Operation. Novye
Tekhnologii, Moscow.
[Lawler, 1979] Lawler, E.L. (1979). Efficient implementation of dynamic programming algorithms for sequencing
problems. CWI Technical report. Stichting Mathematisch Centrum. Mathematische Besliskunde-BW
106/79, 1–16.
[Leon, 1996] Leon, V.J., Peters, B.A. (1996). Replanning and analysis of partial setup strategies in printed circuit
board assembly systems. Int J Flex Manuf Syst 8, 389–411.
[Litl, 1965] Litl, Dzh., Murti, K., Suini, D., Kjerel, K. (1965). Algorithms for Solving the Traveling Salesman
Problem. Economics and Mathematical Methods, vol.1, 94–107 (in Russian).
[Melamed, 1989] Melamed, I.I., Sergeev, S.I., Sigal, I. (1989). The traveling salesman problem. I. Issues in
Theory; II Exact Methods; III Approximate Algorithms. Automation and Remote Control. 50(9), 1147–
1173; 50(10), 1303–1324; 50(11), 1459–1479.
[Schmidt, 1993] Schmidt, G., Strohlein, T. (1993). Relations and Graphs: Discrete Mathematics for Computer
Scientists. EATCS Monographs on Theoretical Computer Science, Springer-Verlag.
[Steiner, 1990] Steiner, G. (1990). On the complexity of dynamic programming for sequencing problems with
precedence constraints. Annals of Operations Research, 26(1), 103–123.
[Tashlykov, 2011] Tashlykov, O.L. (2011). Personnel Dose Costs in the Nuclear Industry. Analysis. Ways to
Decrease. Optimization. LAP LAMBERT Academic Publishing GmbH & Co. RG., Saarbruke.
130